library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.5.1     ✔ tibble    3.2.1
## ✔ lubridate 1.9.4     ✔ tidyr     1.3.1
## ✔ purrr     1.0.2     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(tableone)
library(writexl)
library(gtsummary)
library(survival)
library(ggplot2)
library(readxl)
library(ggpubr)
library(splines)
library(xlsx)
library(lattice)
library(mice)
## 
## Attaching package: 'mice'
## 
## The following object is masked from 'package:stats':
## 
##     filter
## 
## The following objects are masked from 'package:base':
## 
##     cbind, rbind
library(DescTools)
library(survminer)
## 
## Attaching package: 'survminer'
## 
## The following object is masked from 'package:survival':
## 
##     myeloma
library(finalfit)
library(tidycmprsk)
## 
## Attaching package: 'tidycmprsk'
## 
## The following object is masked from 'package:gtsummary':
## 
##     trial
library(sqldf)
## Loading required package: gsubfn
## Loading required package: proto
## Warning in doTryCatch(return(expr), name, parentenv, handler): unable to load shared object '/Library/Frameworks/R.framework/Resources/modules//R_X11.so':
##   dlopen(/Library/Frameworks/R.framework/Resources/modules//R_X11.so, 0x0006): Library not loaded: /opt/X11/lib/libSM.6.dylib
##   Referenced from: <34C5A480-1AC4-30DF-83C9-30A913FC042E> /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/modules/R_X11.so
##   Reason: tried: '/opt/X11/lib/libSM.6.dylib' (no such file), '/System/Volumes/Preboot/Cryptexes/OS/opt/X11/lib/libSM.6.dylib' (no such file), '/opt/X11/lib/libSM.6.dylib' (no such file), '/Library/Frameworks/R.framework/Resources/lib/libSM.6.dylib' (no such file), '/Library/Java/JavaVirtualMachines/jdk-11.0.18+10/Contents/Home/lib/server/libSM.6.dylib' (no such file)
## tcltk DLL is linked to '/opt/X11/lib/libX11.6.dylib'
## Could not load tcltk.  Will use slower R code instead.
## Loading required package: RSQLite
library(ggsurvfit)
library(gridExtra)
## 
## Attaching package: 'gridExtra'
## 
## The following object is masked from 'package:dplyr':
## 
##     combine
library(table1)
## 
## Attaching package: 'table1'
## 
## The following objects are masked from 'package:base':
## 
##     units, units<-
master <- read_xlsx('/Users/to909/Downloads/dataset_updated.xlsx') 
## New names:
## • `` -> `...77`
## • `Cockcroft...80` -> `Cockcroft...81`
## • `study_date...111` -> `study_date...114`
## • `Cockcroft...155` -> `Cockcroft...158`
## • `Cockcroft...268` -> `Cockcroft...271`
## • `study_date...273` -> `study_date...276`
master <- master %>% rename(mGFRadj=`mGFRadj_mL/min/1.73m2`,
                            mGFRunadj=`mGFRunadj_mL/min`,
                            CKDEPIcr21_actualBSA=`CKDEPIcr21_actualBSA_mL/min`,
                            CKDEPIcrcys21_actualBSA=`CKDEPIcrcys21_actualBSA_mL/min`)


# test <- read_xlsx("/Users/to909/Downloads/OncoGFR_MGH_final_0428.xlsx") 

master <- master %>% 
  mutate(bmi_cat=case_when( BMI<18.5 ~ "Underweight",
                            BMI>=18.5 & BMI< 24.9 ~ "Normal Range",
                            BMI>=24.9 & BMI<29.9 ~ "Overweight",
                            BMI>=29.9 ~"Obese")) %>% 
  mutate(bmi_2_cat=ifelse(BMI<25,"<25",">=25"))%>% 
  mutate(sex=gender)


Cockcroft <- function(age,weight_kg,male,cre){
  if(male==0){
    value=(140-age)*weight_kg*0.85/(72*cre)}
  else{
    value=(140-age)*weight_kg/(72*cre)
  }
}

# GET COCKCROFT MEDIAN pre 90 days


master$Cockcroft <-  mapply(Cockcroft,master$age,master$weight,master$male,master$Cr)



master <- master %>% mutate(diff_CKDEPIcr21_measured=mGFRadj-CKDEPIcr21,
                            diff_CKDEPIcysadj_measured=mGFRadj-CKDEPIcysadj,
                            diff_CKDEPIcrcys21_measured=mGFRadj-CKDEPIcrcys21,
                            diff_CKDEPIcysb2mbtpadj_measured=mGFRadj-CKDEPIcysb2mbtpadj,
                            diff_CKDEPI4mark_measured=mGFRadj-CKDEPI4mark,
                            diff_CGadj_measured=mGFRadj-CGadj,
                            diff_Cockcroft_measured=mGFRadj-Cockcroft,
                            diff_CKDEPIb2madj_measured=mGFRadj-CKDEPIb2madj,
                            
                            
                            diff_CKDEPIcrb2madj_measured=mGFRadj-CKDEPIcrb2madj,
                            diff_EKFCcr_measured=mGFRadj-EKFCcr,
                            diff_EKFCcys_measured=mGFRadj-EKFCcys,
                            diff_EKFCcrcys_measured=mGFRadj-EKFCcrcys,

                            diff_CKDEPIbtpadj_measured=mGFRadj-CKDEPIbtpadj,
                            diff_CKDEPIb2mbtpadj_measured=mGFRadj-CKDEPIb2mbtpadj,
                            
                            
                            
                            pct_CKDEPIcysadj_measured=(mGFRadj-CKDEPIcysadj)/CKDEPIcysadj*100,
                            pct_CKDEPIcysadj_measured=(mGFRadj-CKDEPIcysadj)/CKDEPIcysadj*100,
                            pct_CKDEPIcrcys21_measured=(mGFRadj-CKDEPIcrcys21)/CKDEPIcrcys21*100,
                            pct_CKDEPIcysb2mbtpadj_measured=(mGFRadj-CKDEPIcysb2mbtpadj)/CKDEPIcysb2mbtpadj*100,
                            pct_CKDEPI4mark_measured=(mGFRadj-CKDEPI4mark)/CKDEPI4mark*100,
                            pct_CGadj_measured=(mGFRadj-CGadj)/CGadj*100,
                            pct_Cockcroft_measured=(mGFRadj-Cockcroft)/Cockcroft*100,
                            pct_CKDEPIb2madj_measured=(mGFRadj-CKDEPIb2madj)/CKDEPIb2madj*100,
                            
                            pct_CKDEPIcrb2madj_measured=(mGFRadj-CKDEPIcrb2madj)/CKDEPIcrb2madj*100,
                            pct_EKFCcr_measured=(mGFRadj-EKFCcr)/EKFCcr*100,
                            pct_EKFCcys_measured=(mGFRadj-EKFCcys)/EKFCcys*100,
                            pct_EKFCcrcys_measured=(mGFRadj-EKFCcrcys)/EKFCcrcys*100
                            
                            
                            
) 
master <- master %>%
  mutate(race_cat = case_when(
    race == 1 ~ "white",
    race == 2 ~ "black",
    race == 3 ~ "mixed",
    race == 4 ~ "Asian",
    TRUE ~ NA_character_
  ))

master <- master %>%
  mutate(mGFR_cat = case_when(
    mGFRadj >= 90           ~ "≥90",
    mGFRadj >= 60 & mGFRadj <= 89.99 ~ "60-89",
    mGFRadj >= 45 & mGFRadj <= 59.99 ~ "45-59",
    mGFRadj < 45            ~ "<45",
    TRUE ~ NA_character_
  ))

master <- master %>%
  mutate(within_30 = if_else(
    days_between_ct_mGFR >= -30 & days_between_ct_mGFR <= 30,
    1,
    0
  ))

master <- master %>%
  mutate(bias = mGFRadj - CKDEPIcr21)


master$within_30 <- factor(master$within_30, levels = c(0, 1), labels = c("No", "Yes"))
master$sarcopenia_prado <- factor(master$sarcopenia_prado, levels = c(0, 1), labels = c("No Sarcopenia", "Sarcopenia"))
dose_adjustment_level <- function(drug, egfr) {
  if (is.na(egfr)) return(NA)
  # cre_lower <- under dosed
  drug <- tolower(drug)
  
  level <- switch(drug,
    "cisplatin" = {
      if (egfr > 60) 1
      else if (egfr >= 50) 2
      else if (egfr >= 40) 3
      else 5
    },
  
    
    NA  # Return NA if unknown drug
  )
  
  return(level)
}

generate_dose_adjustments <- function(drugs, data, var1, var2) {
  for (drug in drugs) {
    # Calculate levels from var1
    data[[paste0("gold_standard_", drug)]] <- sapply(data[[var1]], function(egfr) {
      dose_adjustment_level(drug, egfr)
    })
    
    # Calculate levels from var2
    data[[paste0("ckdepi_", drug)]] <- sapply(data[[var2]], function(egfr) {
      dose_adjustment_level(drug, egfr)
    })
    
    # Compare and assign label
    data[[paste0("check_", drug)]] <- ifelse(
      data[[paste0("gold_standard_", drug)]] < data[[paste0("ckdepi_", drug)]],
      "under_dosed",
      ifelse(
        data[[paste0("ckdepi_", drug)]] < data[[paste0("gold_standard_", drug)]],
        "over_dosed",
        "same"
      )
    )
  }
  return(data)
}


drugs <- c("cisplatin")


generate_dose_adjustments_new <- function(drugs, data, var1, var2) {
  for (drug in drugs) {
    var1_label <- var1
    var2_label <- var2

    # Calculate levels from var1
    data[[paste0(var1_label, "_", drug)]] <- sapply(data[[var1]], function(egfr) {
      dose_adjustment_level(drug, egfr)
    })

    # Calculate levels from var2
    data[[paste0(var2_label, "_", drug)]] <- sapply(data[[var2]], function(egfr) {
      dose_adjustment_level(drug, egfr)
    })

    # Compare and assign label with updated column name
    data[[paste0("check_", var2_label, "_", drug)]] <- ifelse(
      data[[paste0(var1_label, "_", drug)]] < data[[paste0(var2_label, "_", drug)]],
      "under_dosed",
      ifelse(
        data[[paste0(var2_label, "_", drug)]] < data[[paste0(var1_label, "_", drug)]],
        "over_dosed",
        "same"
      )
    )
  }
  return(data)
}



master_CKDEPIcr21_actualBSA <- generate_dose_adjustments_new(
  drugs = drugs,
  data = master %>% select(protocol,mGFRunadj,CKDEPIcr21_actualBSA),
  var1 = "mGFRunadj",
  var2 = "CKDEPIcr21_actualBSA"
)%>% select(protocol,mGFRunadj_cisplatin,CKDEPIcr21_actualBSA_cisplatin,check_CKDEPIcr21_actualBSA_cisplatin)

master_CKDEPIcrcys21_actualBSA <- generate_dose_adjustments_new(
  drugs = drugs,
  data = master%>% select(protocol,mGFRunadj,CKDEPIcrcys21_actualBSA),
  var1 = "mGFRunadj",
  var2 = "CKDEPIcrcys21_actualBSA"
) %>% select(protocol,CKDEPIcrcys21_actualBSA_cisplatin,check_CKDEPIcrcys21_actualBSA_cisplatin)


dose_check <- master %>% select(protocol,sarcopenia_prado) %>% merge(master_CKDEPIcr21_actualBSA ,all.x = T,by="protocol")%>% merge(master_CKDEPIcrcys21_actualBSA,all.x = T,by="protocol")


names(dose_check)
## [1] "protocol"                               
## [2] "sarcopenia_prado"                       
## [3] "mGFRunadj_cisplatin"                    
## [4] "CKDEPIcr21_actualBSA_cisplatin"         
## [5] "check_CKDEPIcr21_actualBSA_cisplatin"   
## [6] "CKDEPIcrcys21_actualBSA_cisplatin"      
## [7] "check_CKDEPIcrcys21_actualBSA_cisplatin"

Table 1 by sarcopenia groups

table1(~race_cat+mGFR_cat|sarcopenia_prado,caption="By sarcopenia_prado",
       render.continuous = c(.="Mean (SD)", .="Median [Q1,Q3]"),data=master,, footnote="bias = CKDEPIcradj - mGFRadj")
By sarcopenia_prado
No Sarcopenia
(N=308)
Sarcopenia
(N=157)
Overall
(N=465)

bias = CKDEPIcradj - mGFRadj

race_cat
Asian 3 (1.0%) 5 (3.2%) 8 (1.7%)
black 39 (12.7%) 15 (9.6%) 54 (11.6%)
mixed 56 (18.2%) 15 (9.6%) 71 (15.3%)
white 210 (68.2%) 122 (77.7%) 332 (71.4%)
mGFR_cat
<45 13 (4.2%) 22 (14.0%) 35 (7.5%)
≥90 123 (39.9%) 30 (19.1%) 153 (32.9%)
45-59 27 (8.8%) 36 (22.9%) 63 (13.5%)
60-89 145 (47.1%) 69 (43.9%) 214 (46.0%)
test <- master %>% filter(is.na(mGFR_cat)) %>% select(protocol,mGFRadj)

Table cisplatin dose by sarcopenia groups

table1(~factor(CKDEPIcr21_actualBSA_cisplatin)+check_CKDEPIcr21_actualBSA_cisplatin+factor(CKDEPIcrcys21_actualBSA_cisplatin)+check_CKDEPIcrcys21_actualBSA_cisplatin
       |sarcopenia_prado,caption="By sarcopenia_prado",footnote="cisplatin = 
      (egfr > 60) 1
      (egfr >= 50) 2
      (egfr >= 40) 3
      else 5
    ",test = TRUE,
       render.continuous = c(.="Mean (SD)", .="Median [Q1,Q3]"),data=dose_check)
By sarcopenia_prado
No Sarcopenia
(N=308)
Sarcopenia
(N=157)
Overall
(N=465)

cisplatin = (egfr > 60) 1 (egfr >= 50) 2 (egfr >= 40) 3 else 5

factor(CKDEPIcr21_actualBSA_cisplatin)
1 286 (92.9%) 132 (84.1%) 418 (89.9%)
2 8 (2.6%) 6 (3.8%) 14 (3.0%)
3 9 (2.9%) 10 (6.4%) 19 (4.1%)
5 5 (1.6%) 9 (5.7%) 14 (3.0%)
check_CKDEPIcr21_actualBSA_cisplatin
over_dosed 27 (8.8%) 40 (25.5%) 67 (14.4%)
same 278 (90.3%) 115 (73.2%) 393 (84.5%)
under_dosed 3 (1.0%) 2 (1.3%) 5 (1.1%)
factor(CKDEPIcrcys21_actualBSA_cisplatin)
1 273 (88.6%) 112 (71.3%) 385 (82.8%)
2 16 (5.2%) 19 (12.1%) 35 (7.5%)
3 11 (3.6%) 13 (8.3%) 24 (5.2%)
5 8 (2.6%) 13 (8.3%) 21 (4.5%)
check_CKDEPIcrcys21_actualBSA_cisplatin
over_dosed 21 (6.8%) 26 (16.6%) 47 (10.1%)
same 273 (88.6%) 120 (76.4%) 393 (84.5%)
under_dosed 14 (4.5%) 11 (7.0%) 25 (5.4%)

Table within 30 days only

table1(~bias|within_30,caption="By within_30",
       render.continuous = c(.="Mean (SD)", .="Median [Q1,Q3]"),data=master,, footnote="bias = CKDEPIcradj - mGFRadj")
By within_30
No
(N=127)
Yes
(N=338)
Overall
(N=465)

bias = CKDEPIcradj - mGFRadj

bias
Mean (SD) -8.96 (13.1) -12.1 (15.0) -11.2 (14.6)
Median [Q1,Q3] -10.1 [-16.3,0.695] -12.0 [-21.9,-2.25] -11.5 [-20.7,-0.970]
wilcox.test(bias ~ within_30, data = master)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  bias by within_30
## W = 24238, p-value = 0.03164
## alternative hypothesis: true location shift is not equal to 0

Table within_30 by sarcopenia groups

master <- master %>% mutate(within_30_Sarcopenia=case_when(within_30=="Yes" & sarcopenia_prado=="Sarcopenia" ~ "<30 days Sarcopenia",
                                                 within_30=="No" & sarcopenia_prado=="Sarcopenia" ~ ">=30 days Sarcopenia",
                                                 within_30=="Yes" & sarcopenia_prado=="No Sarcopenia" ~ "<30 days No Sarcopenia",
                                                 within_30=="No" & sarcopenia_prado=="No Sarcopenia" ~ ">=30 days No Sarcopenia",
                                                 ))

table1(~bias|within_30_Sarcopenia,caption="By within_30",
       render.continuous = c(.="Mean (SD)", .="Median [Q1,Q3]"),data=master, footnote="bias = CKDEPIcradj - mGFRadj")
By within_30
<30 days No Sarcopenia
(N=227)
<30 days Sarcopenia
(N=111)
>=30 days No Sarcopenia
(N=81)
>=30 days Sarcopenia
(N=46)
Overall
(N=465)

bias = CKDEPIcradj - mGFRadj

bias
Mean (SD) -9.25 (14.4) -17.9 (14.7) -6.57 (13.2) -13.2 (11.9) -11.2 (14.6)
Median [Q1,Q3] -9.07 [-18.9,0.330] -17.6 [-26.7,-7.00] -5.58 [-14.8,2.62] -13.2 [-18.6,-5.96] -11.5 [-20.7,-0.970]

Bland-Altman Plot 1 mGFRadj - CKDEPIcr21

library(ggplot2)

# Calculate difference and average
master$average <- (master$mGFRadj + master$CKDEPIcr21) / 2
master$difference <- master$mGFRadj - master$CKDEPIcr21

# Calculate mean difference and limits of agreement
mean_diff <- mean(master$difference, na.rm = TRUE)
loa <- sd(master$difference, na.rm = TRUE) * 1.96
upper_loa <- mean_diff + loa
lower_loa <- mean_diff - loa

# Create Bland-Altman plot
ggplot(master, aes(x = average, y = difference)) +
  geom_point(aes(color = sarcopenia_prado), alpha = 0.6) +
  scale_color_manual(values = c("Sarcopenia" = "#FF0000", "No Sarcopenia" = "#007FFF"),name = "Sarcopenia group") +
  geom_hline(yintercept = 0, color = "blue", linetype = "solid") +
  geom_hline(yintercept = c(-15, 15), color = "red", linetype = "dashed") +
  geom_hline(yintercept = c(-30, -60), color = "black", linetype = "dotted") +
  geom_hline(yintercept = c(30, 60), color = "black", linetype = "dotted") +
  scale_y_continuous(
    limits = c(-60, 60),
    breaks = c(-60, -30, -15, 0, 15, 30, 60)
  ) +
  labs(
    title = "Bland-Altman Plot",
    x = bquote("Average GFR (mL/min per 1.73 m"^2*")"),,
    y = "Measured GFR minus eGFR"
  )+
  theme_minimal()+scale_x_continuous(limits = c(0, 150), breaks = seq(0, 150, by = 50))

Bland-Altman Plot 2 mGFRadj - CKDEPIcysadj

library(ggplot2)

# Calculate difference and average
master$average <- (master$mGFRadj + master$CKDEPIcysadj) / 2
master$difference <- master$mGFRadj - master$CKDEPIcysadj

# Calculate mean difference and limits of agreement
mean_diff <- mean(master$difference, na.rm = TRUE)
loa <- sd(master$difference, na.rm = TRUE) * 1.96
upper_loa <- mean_diff + loa
lower_loa <- mean_diff - loa

# Create Bland-Altman plot
ggplot(master, aes(x = average, y = difference)) +
  geom_point(aes(color = sarcopenia_prado), alpha = 0.6) +
  scale_color_manual(values = c("Sarcopenia" = "#FF0000", "No Sarcopenia" = "#007FFF"),name = "Sarcopenia group") +
  geom_hline(yintercept = 0, color = "blue", linetype = "solid") +
  geom_hline(yintercept = c(-15, 15), color = "red", linetype = "dashed") +
  geom_hline(yintercept = c(-30, -60), color = "black", linetype = "dotted") +
  geom_hline(yintercept = c(30, 60), color = "black", linetype = "dotted") +
  scale_y_continuous(
    limits = c(-60, 60),
    breaks = c(-60, -30, -15, 0, 15, 30, 60)
  ) +
  labs(
    title = "Bland-Altman Plot",
    x = bquote("Average GFR (mL/min per 1.73 m"^2*")"),,
    y = "Measured GFR minus eGFR"
  ) +
  theme_minimal()+scale_x_continuous(limits = c(0, 150), breaks = seq(0, 150, by = 50))

Bland-Altman Plot 3 mGFRadj - CKDEPIcrcys21

library(ggplot2)

# Calculate difference and average
master$average <- (master$mGFRadj + master$CKDEPIcrcys21) / 2
master$difference <- master$mGFRadj - master$CKDEPIcrcys21

# Calculate mean difference and limits of agreement
mean_diff <- mean(master$difference, na.rm = TRUE)
loa <- sd(master$difference, na.rm = TRUE) * 1.96
upper_loa <- mean_diff + loa
lower_loa <- mean_diff - loa

# Create Bland-Altman plot
ggplot(master, aes(x = average, y = difference)) +
  geom_point(aes(color = sarcopenia_prado), alpha = 0.6) +
  scale_color_manual(values = c("Sarcopenia" = "#FF0000", "No Sarcopenia" = "#007FFF"),name = "Sarcopenia group") +
  geom_hline(yintercept = 0, color = "blue", linetype = "solid") +
  geom_hline(yintercept = c(-15, 15), color = "red", linetype = "dashed") +
  geom_hline(yintercept = c(-30, -60), color = "black", linetype = "dotted") +
  geom_hline(yintercept = c(30, 60), color = "black", linetype = "dotted") +
  scale_y_continuous(
    limits = c(-60, 60),
    breaks = c(-60, -30, -15, 0, 15, 30, 60)
  ) +
  labs(
    title = "Bland-Altman Plot",
    x = bquote("Average GFR (mL/min per 1.73 m"^2*")"),,
    y = "Measured GFR minus eGFR"
  ) +
  theme_minimal()+scale_x_continuous(limits = c(0, 150), breaks = seq(0, 150, by = 50))

Bland-Altman Plot 4 mGFRadj - CKDEPIcrb2madj

library(ggplot2)

# Calculate difference and average
master$average <- (master$mGFRadj + master$CKDEPIcrb2madj) / 2
master$difference <- master$mGFRadj - master$CKDEPIcrb2madj

# Calculate mean difference and limits of agreement
mean_diff <- mean(master$difference, na.rm = TRUE)
loa <- sd(master$difference, na.rm = TRUE) * 1.96
upper_loa <- mean_diff + loa
lower_loa <- mean_diff - loa

# Create Bland-Altman plot
ggplot(master, aes(x = average, y = difference)) +
  geom_point(aes(color = sarcopenia_prado), alpha = 0.6) +
  scale_color_manual(values = c("Sarcopenia" = "#FF0000", "No Sarcopenia" = "#007FFF"),name = "Sarcopenia group") +
  geom_hline(yintercept = 0, color = "blue", linetype = "solid") +
  geom_hline(yintercept = c(-15, 15), color = "red", linetype = "dashed") +
  geom_hline(yintercept = c(-30, -60), color = "black", linetype = "dotted") +
  geom_hline(yintercept = c(30, 60), color = "black", linetype = "dotted") +
  scale_y_continuous(
    limits = c(-60, 60),
    breaks = c(-60, -30, -15, 0, 15, 30, 60)
  ) +
  labs(
    title = "Bland-Altman Plot",
    x = bquote("Average GFR (mL/min per 1.73 m"^2*")"),,
    y = "Measured GFR minus eGFR"
  ) +
  theme_minimal()+scale_x_continuous(limits = c(0, 150), breaks = seq(0, 150, by = 50))

Bland-Altman Plot 5 mGFRadj - CKDEPIcysb2mbtpadj

library(ggplot2)

# Calculate difference and average
master$average <- (master$mGFRadj + master$CKDEPIcysb2mbtpadj) / 2
master$difference <- master$mGFRadj - master$CKDEPIcysb2mbtpadj

# Calculate mean difference and limits of agreement
mean_diff <- mean(master$difference, na.rm = TRUE)
loa <- sd(master$difference, na.rm = TRUE) * 1.96
upper_loa <- mean_diff + loa
lower_loa <- mean_diff - loa

# Create Bland-Altman plot
ggplot(master, aes(x = average, y = difference)) +
  geom_point(aes(color = sarcopenia_prado), alpha = 0.6) +
  scale_color_manual(values = c("Sarcopenia" = "#FF0000", "No Sarcopenia" = "#007FFF"),name = "Sarcopenia group") +
  geom_hline(yintercept = 0, color = "blue", linetype = "solid") +
  geom_hline(yintercept = c(-15, 15), color = "red", linetype = "dashed") +
  geom_hline(yintercept = c(-30, -60), color = "black", linetype = "dotted") +
  geom_hline(yintercept = c(30, 60), color = "black", linetype = "dotted") +
  scale_y_continuous(
    limits = c(-60, 60),
    breaks = c(-60, -30, -15, 0, 15, 30, 60)
  ) +
  labs(
    title = "Bland-Altman Plot",
    x = bquote("Average GFR (mL/min per 1.73 m"^2*")"),,
    y = "Measured GFR minus eGFR"
  ) +
  theme_minimal()+scale_x_continuous(limits = c(0, 150), breaks = seq(0, 150, by = 50))

Bland-Altman Plot 6 mGFRadj - CKDEPI4mark

library(ggplot2)

# Calculate difference and average
master$average <- (master$mGFRadj + master$CKDEPI4mark) / 2
master$difference <- master$mGFRadj - master$CKDEPI4mark

# Calculate mean difference and limits of agreement
mean_diff <- mean(master$difference, na.rm = TRUE)
loa <- sd(master$difference, na.rm = TRUE) * 1.96
upper_loa <- mean_diff + loa
lower_loa <- mean_diff - loa


# Create Bland-Altman plot
ggplot(master, aes(x = average, y = difference)) +
  geom_point(aes(color = sarcopenia_prado), alpha = 0.6) +
  scale_color_manual(values = c("Sarcopenia" = "#FF0000", "No Sarcopenia" = "#007FFF"),name = "Sarcopenia group") +
  geom_hline(yintercept = 0, color = "blue", linetype = "solid") +
  geom_hline(yintercept = c(-15, 15), color = "red", linetype = "dashed") +
  geom_hline(yintercept = c(-30, -60), color = "black", linetype = "dotted") +
  geom_hline(yintercept = c(30, 60), color = "black", linetype = "dotted") +
  scale_y_continuous(
    limits = c(-60, 60),
    breaks = c(-60, -30, -15, 0, 15, 30, 60)
  ) +
  labs(
    title = "Bland-Altman Plot",
    x = bquote("Average GFR (mL/min per 1.73 m"^2*")"),,
    y = "Measured GFR minus eGFR"
  ) +
  theme_minimal()+scale_x_continuous(limits = c(0, 150), breaks = seq(0, 150, by = 50))