Toward the standardization of big datasets of urine output for AKI analysis: A multicenter validation study

 

This is a full reproduction for the results of the validation cohort (AUMCdb).

 

Citation

Hasidim, A.A., Klein, M.A., Ben Shitrit, I. et al. Toward the standardization of big datasets of urine output for AKI analysis: a multicenter validation study. Sci Rep 15, 20009 (2025). https://doi.org/10.1038/s41598-025-95535-4

---

Population and Study’s Sample

Total I/M-CU stays in aumc:

print("Medium and intensive care:")

## [1] "Medium and intensive care:"

dbGetQuery(con, "SELECT count(*) FROM `original.admissions`")

## # A tibble: 1 × 1
##     f0_
##   <int>
## 1 23106

print("Intensive care only:")

## [1] "Intensive care only:"

dbGetQuery(con, "SELECT count(*) FROM `original.admissions` WHERE lower(location) != 'mc'")

## # A tibble: 1 × 1
##     f0_
##   <int>
## 1 18386

ICU stays with UO records (eligible):

print("Raw UO per service type:")

## [1] "Raw UO per service type:"

table(raw_uo$SERVICE)

## 
##      IC   IC&MC      MC   MC&IC 
## 1189758  374109  106537    9666

print("Total eligible ICU admissions:")

## [1] "Total eligible ICU admissions:"

raw_uo %>% filter(SERVICE != "MC") %>%
  {n_distinct(.$STAY_ID)}

## [1] 18147

Patients with UO records:

raw_uo %>% filter(SERVICE != "MC") %>%
  {n_distinct(.$SUBJECT_ID)}

## [1] 16344

 

Count all UO records (before exclusions):

all_rows_count <- nrow(raw_uo %>% filter(SERVICE != "MC"))
all_rows_count

## [1] 1573533

 

Figure for Raw UO Records Data

Frequency of urine output charting by source

bS3a <- raw_uo %>% 
  mutate(
    LABEL = case_when(
      LABEL == "GU Irrigant/Urine Volume Out" ~ "GU Irrig. Out",
      LABEL == "GU Irrigant Volume In" ~ "GU Irrig. In",
      .default = LABEL
    )
  ) %>% count(LABEL, sort = TRUE) %>%
   ggplot(aes(x=reorder(LABEL, -n), y=n)) +
   geom_bar(stat="identity") +
   xlab("") +
   ylab("") +
   geom_bar(stat="identity", fill="steelblue") +
   geom_text(aes(label=n), vjust=-0.6, color="black", size=3) +
   theme_classic() +
   theme(axis.text.y=element_blank()) +
   theme(axis.text.x = element_text(angle = 45, hjust = 1))

bS3a

 

Check for Duplications

Distinctiveness

First, we are basing distinctive rows in the raw UO data.

Count distinct raw rows:

distinct_time_item_patient_rows_count <- raw_uo %>% 
  select(-VALUE, 
         -SERVICE, 
         -LABEL) %>% 
  n_distinct()

distinct_time_item_patient_rows_count

## [1] 1680070

Conclusion: the original raw query does not have duplicates (all rows are distinct by all columns)

 

Simultaneous Charting

raw_uo_excluions_duplicates$same_value <- as.factor(raw_uo_excluions_duplicates$same_value)
raw_uo_excluions_duplicates$label <- as.factor(raw_uo_excluions_duplicates$label)
raw_uo_excluions_duplicates$label <- factor(raw_uo_excluions_duplicates$label, levels = as.factor(names(sort(table(raw_uo_excluions_duplicates$label),
                                  decreasing = TRUE))))

S1_a <- raw_uo_excluions_duplicates %>%
  select(same_value, label) %>%
  tbl_summary(by = same_value)

S1_a

Characteristic	Different volume N = 6,1801	Equal volume N = 5731
label
Foley,R Nephrostomy	1,200 (19%)	189 (33%)
Foley,L Nephrostomy	967 (16%)	102 (18%)
Foley,L Nephrostomy,R Nephrostomy	942 (15%)	40 (7.0%)
Foley,Suprapubic	701 (11%)	37 (6.5%)
L Nephrostomy,R Nephrostomy	421 (6.8%)	50 (8.7%)
Suprapubic,R Nephrostomy	331 (5.4%)	28 (4.9%)
Foley,Suprapubic,L Nephrostomy,R Nephrostomy,L Ureteral Stent,R Ureteral Stent	339 (5.5%)	0 (0%)
Foley,Void	182 (2.9%)	10 (1.7%)
Foley,Ileoconduit	159 (2.6%)	21 (3.7%)
Foley,Incontinence (urine leakage)	159 (2.6%)	17 (3.0%)
Suprapubic,Ileoconduit	131 (2.1%)	18 (3.1%)
Ileoconduit,L Nephrostomy,R Nephrostomy	93 (1.5%)	2 (0.3%)
Foley,Suprapubic,L Nephrostomy,R Nephrostomy	61 (1.0%)	3 (0.5%)
Foley,Suprapubic,Ileoconduit	61 (1.0%)	0 (0%)
Suprapubic,Void	59 (1.0%)	1 (0.2%)
Foley,L Ureteral Stent,R Ureteral Stent	57 (0.9%)	1 (0.2%)
Void,Incontinence (urine leakage)	38 (0.6%)	7 (1.2%)
Void,L Nephrostomy	35 (0.6%)	4 (0.7%)
Ileoconduit,R Ureteral Stent	24 (0.4%)	11 (1.9%)
Ileoconduit,R Nephrostomy	32 (0.5%)	0 (0%)
Incontinence (urine leakage),Ileoconduit	29 (0.5%)	1 (0.2%)
Foley,L Ureteral Stent	25 (0.4%)	3 (0.5%)
L Ureteral Stent,R Ureteral Stent	7 (0.1%)	15 (2.6%)
Ileoconduit,L Nephrostomy,R Nephrostomy,R Ureteral Stent	16 (0.3%)	0 (0%)
Suprapubic,Incontinence (urine leakage)	12 (0.2%)	3 (0.5%)
Suprapubic,R Ureteral Stent	12 (0.2%)	1 (0.2%)
Void,R Nephrostomy	11 (0.2%)	1 (0.2%)
L Nephrostomy,R Nephrostomy,R Ureteral Stent	10 (0.2%)	0 (0%)
Suprapubic,L Nephrostomy	7 (0.1%)	2 (0.3%)
Foley,Ileoconduit,L Nephrostomy	8 (0.1%)	0 (0%)
Foley,L Nephrostomy,R Nephrostomy,L Ureteral Stent,R Ureteral Stent	6 (<0.1%)	0 (0%)
Foley,Suprapubic,L Nephrostomy	5 (<0.1%)	0 (0%)
Ileoconduit,L Nephrostomy	4 (<0.1%)	1 (0.2%)
Incontinence (urine leakage),L Nephrostomy	5 (<0.1%)	0 (0%)
L Nephrostomy,R Nephrostomy,L Ureteral Stent,R Ureteral Stent	5 (<0.1%)	0 (0%)
Suprapubic,L Nephrostomy,R Nephrostomy	4 (<0.1%)	1 (0.2%)
Ileoconduit,R Nephrostomy,R Ureteral Stent	4 (<0.1%)	0 (0%)
Suprapubic,L Nephrostomy,R Nephrostomy,L Ureteral Stent,R Ureteral Stent	3 (<0.1%)	0 (0%)
Void,Ileoconduit	2 (<0.1%)	1 (0.2%)
Foley,Suprapubic,Incontinence (urine leakage)	2 (<0.1%)	0 (0%)
Foley,Suprapubic,R Nephrostomy,L Ureteral Stent,R Ureteral Stent	2 (<0.1%)	0 (0%)
Foley,Void,Incontinence (urine leakage)	2 (<0.1%)	0 (0%)
Foley,Suprapubic,L Nephrostomy,R Nephrostomy,L Ureteral Stent	1 (<0.1%)	0 (0%)
Foley,Suprapubic,L Nephrostomy,R Nephrostomy,R Ureteral Stent	1 (<0.1%)	0 (0%)
Foley,Suprapubic,R Nephrostomy,R Ureteral Stent	1 (<0.1%)	0 (0%)
Foley,Suprapubic,Void	1 (<0.1%)	0 (0%)
Foley,Void,Ileoconduit	0 (0%)	1 (0.2%)
Ileoconduit,L Nephrostomy,R Ureteral Stent	1 (<0.1%)	0 (0%)
R Nephrostomy,R Ureteral Stent	0 (0%)	1 (0.2%)
Suprapubic,Incontinence (urine leakage),Ileoconduit	0 (0%)	1 (0.2%)
Suprapubic,Void,Incontinence (urine leakage)	1 (<0.1%)	0 (0%)
Void,L Nephrostomy,R Nephrostomy	1 (<0.1%)	0 (0%)
1 n (%)

      Show full SQL query —–>

WITH aa AS (
        SELECT
            STAY_ID,
            CHARTTIME
        FROM
            `aumc_uo_and_aki.a_urine_output_raw` uo
        GROUP BY
            STAY_ID,
            CHARTTIME
    ),
    ab AS (
        SELECT
            a.*,
            b.ITEMID,
            b.VALUE,
            b.LABEL
        FROM
            aa a
            LEFT JOIN `aumc_uo_and_aki.a_urine_output_raw` b ON b.STAY_ID = a.STAY_ID
            AND b.CHARTTIME = a.CHARTTIME
        ORDER BY
            STAY_ID,
            CHARTTIME,
            ITEMID
    ),
    ac AS (
        SELECT
            STAY_ID,
            CHARTTIME,
            STRING_AGG(label) label,
            COUNT(STAY_ID) COUNT,
            IF(
                MIN(VALUE) = MAX(VALUE),
                "Equal volume",
                "Different volume"
            ) AS same_value
        FROM
            ab
        GROUP BY
            STAY_ID,
            CHARTTIME
    ),
    ad AS (
        SELECT
            COUNT(CHARTTIME) COUNT,
            label,
            same_value
        FROM
            ac
        WHERE
            COUNT > 1
        GROUP BY
            label,
            same_value
    )
SELECT
    *
FROM
    ac
WHERE
    COUNT > 1

In conclusion, most of the records have different values, and thus human error in duplicate record-keeping is not likely.

 

Exclusion

ICU type exclusion:

dbGetQuery(con, statement = read_file('sql/service_type_exclusion.sql'))

## # A tibble: 1 × 2
##   icu_stays UO_records
##       <int>      <int>
## 1      4641     106537

      Show full SQL query —–>

SELECT
    COUNT(DISTINCT STAY_ID) icu_stays,
    COUNT(STAY_ID) UO_records
FROM
    `aumc_uo_and_aki.a_urine_output_raw`
WHERE
    lower(SERVICE) = "mc"

 

Uretral stent exclusion:

dbGetQuery(con, statement = read_file('sql/ure_stent_exclusion.sql'))

## # A tibble: 1 × 2
##   icu_stays UO_records
##       <int>      <int>
## 1        10       3804

      Show full SQL query —–>

SELECT
    COUNT(DISTINCT STAY_ID) icu_stays,
    COUNT(STAY_ID) UO_records
FROM
    `aumc_uo_and_aki.a_urine_output_raw`
WHERE
    STAY_ID IN (
        SELECT
            admissionid AS STAY_ID
        FROM
            `original.numericitems`
        WHERE
            ITEMID IN (19921, 19922) -- Urethral stent
        GROUP BY
            admissionid
    )
    AND lower(SERVICE) != "mc"

 

Urine incontinence exclusion:

dbGetQuery(con, statement = read_file('sql/Incontinence_exclusion.sql'))

## # A tibble: 1 × 2
##   icu_stays UO_records
##       <int>      <int>
## 1       250      62136

      Show full SQL query —–>

SELECT
    COUNT(DISTINCT STAY_ID) icu_stays,
    COUNT(STAY_ID) UO_records
FROM
    `aumc_uo_and_aki.a_urine_output_raw`
WHERE
    STAY_ID IN (
        SELECT
            admissionid AS STAY_ID
        FROM
            `original.numericitems`
        WHERE
            ITEMID IN (8800) -- Incontinence (Urine leakage)
        GROUP BY
            admissionid
    )
    AND STAY_ID NOT IN (
        SELECT
            admissionid AS STAY_ID
        FROM
            `original.numericitems`
        WHERE
            ITEMID IN (19921, 19922) -- Urethral stent
        GROUP BY
            admissionid
    )
    AND lower(SERVICE) != "mc"

 

Not passing UO sanity check:

dbGetQuery(con, statement = read_file('sql/sanity.sql'))

## # A tibble: 1 × 1
##   UO_records
##        <int>
## 1          4

      Show full SQL query —–>

SELECT
  COUNT(STAY_ID) UO_records
FROM
  `aumc_uo_and_aki.a_urine_output_raw`
WHERE
  STAY_ID NOT IN (
    SELECT
      admissionid AS STAY_ID
    FROM
      `original.numericitems`
    WHERE
      ITEMID IN (8800) -- Incontinence (Urine leakage)
    GROUP BY
      admissionid
  )
  AND STAY_ID NOT IN (
    SELECT
      admissionid AS STAY_ID
    FROM
      `original.numericitems`
    WHERE
      ITEMID IN (19921, 19922) -- Urethral stent
    GROUP BY
      admissionid
  )
  AND LOWER(SERVICE) != "mc"
  AND (
    VALUE > 5000
    OR VALUE < 0
  )

 

Total raw urine output after exclusion (“eligible records, before collection time”):

nrow(raw_uo_eligible)

## [1] 1507589

      Show full SQL query —–>

SELECT
  *
FROM
  `aumc_uo_and_aki.a_urine_output_raw`
WHERE
  STAY_ID NOT IN (
    SELECT
      admissionid AS STAY_ID
    FROM
      `original.numericitems`
    WHERE
      ITEMID IN (19922, 19921, 8800) -- Urethral stent, Urine Incontinence
    GROUP BY
      admissionid
  )
  AND NOT (
    VALUE > 5000
    OR VALUE < 0
  )
  AND lower(SERVICE) != "mc"

 

Total icu stays after exclusion (“eligible icu admissions, before collection time”):

n_distinct(raw_uo_eligible$STAY_ID)

## [1] 17887

Exclusion of first volume in each compartment per ICU stay:

uo_rate_including_null_collection_period %>%
  filter(STAY_ID %in% raw_uo_eligible$STAY_ID) %>%
  filter(is.na(TIME_INTERVAL)) %>%
  nrow()

## [1] 18062

 

UO records with time intervals:

uo_rate_including_null_collection_period %>%
  filter(STAY_ID %in% raw_uo_eligible$STAY_ID) %>%
  drop_na(TIME_INTERVAL) %>%
  nrow()

## [1] 1489527

Count UO records by anatomical compartment:

uo_rate %>% 
  mutate(agg_group = case_when(SOURCE == "Foley" |
                                 SOURCE == "Condom Cath" |
                                 SOURCE == "Straight Cath" |
                                 SOURCE == "Suprapubic" |
                                 SOURCE == "Void" ~ "Urinary bladder",
                               TRUE ~ SOURCE)
  ) %>%
           group_by(agg_group) %>%
   dplyr::summarise(N = n()
  )

## # A tibble: 4 × 2
##   agg_group             N
##   <chr>             <int>
## 1 Ileoconduit        4749
## 2 L Nephrostomy      2404
## 3 R Nephrostomy      3423
## 4 Urinary bladder 1478951

ICU stays with UO records with time intervals:

print("all ICU stays:")

## [1] "all ICU stays:"

uo_rate_including_null_collection_period %>%
  {n_distinct(.$STAY_ID)}

## [1] 17887

print("ICU stays with time intervals:")

## [1] "ICU stays with time intervals:"

uo_rate_including_null_collection_period %>%
                     drop_na(TIME_INTERVAL) %>%
  {n_distinct(.$STAY_ID)}

## [1] 17736

print("ICU stays with UO records that does not  have time interval (no previous UO record in the same compartment:")

## [1] "ICU stays with UO records that does not  have time interval (no previous UO record in the same compartment:"

uo_rate_including_null_collection_period %>%
                     filter(is.na(TIME_INTERVAL)) %>%
  {n_distinct(.$STAY_ID)}

## [1] 17887

print("ICU stays dropped due to no UO records with time intervalst (no previous UO record in the same compartment:")

## [1] "ICU stays dropped due to no UO records with time intervalst (no previous UO record in the same compartment:"

(uo_rate_including_null_collection_period %>%
                     filter(is.na(TIME_INTERVAL))) %>%
  filter(!(STAY_ID %in% (uo_rate_including_null_collection_period %>%
                     drop_na(TIME_INTERVAL))$STAY_ID)) %>%
  {n_distinct(.$STAY_ID)}

## [1] 151

Count total ICU days of UO monitoring for icu stays with time intervals:

hourly_uo %>%
  filter(STAY_ID %in%
           (uo_rate_including_null_collection_period %>%
                     drop_na(TIME_INTERVAL))$STAY_ID) %>%
  nrow() / 24

## [1] 89959.58

Hours of UO monitoring:

print("Eligible hours of UO monitoring:")

## [1] "Eligible hours of UO monitoring:"

hourly_uo %>%
  filter(STAY_ID %in%
           (
             uo_rate_including_null_collection_period %>%
               drop_na(TIME_INTERVAL)
           )$STAY_ID) %>%
  nrow()

## [1] 2159030

print("Valid hourly-adjusted UO monitoring hours:")

## [1] "Valid hourly-adjusted UO monitoring hours:"

hourly_uo %>%
  filter(STAY_ID %in%
           (
             uo_rate_including_null_collection_period %>%
               drop_na(TIME_INTERVAL)
           )$STAY_ID) %>%
  drop_na(HOURLY_WEIGHTED_MEAN_RATE) %>%
  nrow()

## [1] 2140972

print("ICU stays with valid hourly-adjusted UO monitoring hours:")

## [1] "ICU stays with valid hourly-adjusted UO monitoring hours:"

hourly_uo %>%
  filter(STAY_ID %in%
           (
             uo_rate_including_null_collection_period %>%
               drop_na(TIME_INTERVAL)
           )$STAY_ID) %>%
  drop_na(HOURLY_WEIGHTED_MEAN_RATE) %>%
  {n_distinct(.$STAY_ID)}

## [1] 17728

Proportion of valid hours covered out of included hours of uo monitoring:

hourly_uo %>%
  filter(STAY_ID %in% raw_uo_eligible$STAY_ID) %>%
  drop_na(HOURLY_WEIGHTED_MEAN_RATE) %>%
  nrow() / hourly_uo %>%
  filter(STAY_ID %in% raw_uo_eligible$STAY_ID) %>%
  nrow()

## [1] 0.9915644

 

Hourly-adjusted UO with admission weights:

print("ICU stays with weight at admission and valid hourly-adjusted UO:")

## [1] "ICU stays with weight at admission and valid hourly-adjusted UO:"

hourly_uo %>%
  filter(STAY_ID %in% raw_uo_eligible$STAY_ID) %>%
  drop_na(HOURLY_WEIGHTED_MEAN_RATE) %>%
  drop_na(WEIGHT_ADMIT) %>%
  {n_distinct(.$STAY_ID)}

## [1] 17136

print("Valid hourly-adjusted UO monitoring hours with weight at admission:")

## [1] "Valid hourly-adjusted UO monitoring hours with weight at admission:"

hourly_uo %>%
  filter(STAY_ID %in% raw_uo_eligible$STAY_ID) %>%
  drop_na(HOURLY_WEIGHTED_MEAN_RATE) %>%
  drop_na(WEIGHT_ADMIT) %>%
  nrow()

## [1] 2063720

print("ICU stays with valid weight (25 <= kg <=300) and valid hourly-adjusted UO:")

## [1] "ICU stays with valid weight (25 <= kg <=300) and valid hourly-adjusted UO:"

hourly_uo %>%
  filter(STAY_ID %in% raw_uo_eligible$STAY_ID) %>%
  drop_na(HOURLY_WEIGHTED_MEAN_RATE) %>%
  drop_na(WEIGHT_ADMIT) %>%
  filter(WEIGHT_ADMIT <= 300,
         WEIGHT_ADMIT >= 25) %>%
  {n_distinct(.$STAY_ID)}

## [1] 17136

print("Valid hourly-adjusted UO monitoring hours with valid weights:")

## [1] "Valid hourly-adjusted UO monitoring hours with valid weights:"

hourly_uo %>%
  filter(STAY_ID %in% raw_uo_eligible$STAY_ID) %>%
  drop_na(HOURLY_WEIGHTED_MEAN_RATE) %>%
  drop_na(WEIGHT_ADMIT) %>%
  filter(WEIGHT_ADMIT <= 300,
         WEIGHT_ADMIT >= 25) %>%
  nrow()

## [1] 2063720

ICU stays with calculated KDIGO-UO staging (at least six consecutive hours with valid charting of hourly-adjusted UO)

print("number of ICU stays with valid KDIGO-UO staging")

## [1] "number of ICU stays with valid KDIGO-UO staging"

kdigo_uo_stage %>% 
  filter(stay_id %in% raw_uo_eligible$STAY_ID) %>%
  filter(WEIGHT_ADMIT <= 300,
         WEIGHT_ADMIT >= 25) %>%
  {n_distinct(.$stay_id)}

## [1] 16869

print("number of included hours with valid KDIGO-UO staging")

## [1] "number of included hours with valid KDIGO-UO staging"

kdigo_uo_stage %>% 
  filter(stay_id %in% raw_uo_eligible$STAY_ID) %>%
  filter(WEIGHT_ADMIT <= 300,
         WEIGHT_ADMIT >= 25) %>%
  nrow()

## [1] 1978407

eligible first ICU admission to each patient for AKI analysis

first_icu_stay <-
  dbGetQuery(
    con,
    "SELECT patientid, admissionid FIRST_STAY_ID_IN_PATIENT FROM `original.admissions`
where admissioncount = 1"
  )

print("number of first ICU stays with valid KDIGO-UO staging")

## [1] "number of first ICU stays with valid KDIGO-UO staging"

kdigo_uo_stage %>% 
  filter(stay_id %in% first_icu_stay$FIRST_STAY_ID_IN_PATIENT
         & stay_id %in% raw_uo_eligible$STAY_ID) %>%
  filter(WEIGHT_ADMIT <= 300,
         WEIGHT_ADMIT >= 25) %>%
  {n_distinct(.$stay_id)}

## [1] 14932

print("number of valid hourly KDIGO-UO staging in first ICU stays")

## [1] "number of valid hourly KDIGO-UO staging in first ICU stays"

kdigo_uo_stage %>% 
  filter(stay_id %in% first_icu_stay$FIRST_STAY_ID_IN_PATIENT
         & stay_id %in% raw_uo_eligible$STAY_ID) %>%
  filter(WEIGHT_ADMIT <= 300,
         WEIGHT_ADMIT >= 25) %>%
  nrow()

## [1] 1558176

included ICU admissions for AKI analysis

print("number of first ICU stays (after exclusion criteria) with valid KDIGO-UO staging for the first 24-hours ins ICU stay")

## [1] "number of first ICU stays (after exclusion criteria) with valid KDIGO-UO staging for the first 24-hours ins ICU stay"

aki_epi <- akis_all_long %>%
  filter(group == "newcons") %>%
  drop_na(prevalnce_admit) %>%
  transmute(STAY_ID,
           first_kdigo_uo = first_stage,
         max_uo_stage = max_stage)
  
kdigo_uo_stage %>% 
  filter(stay_id %in% aki_epi$STAY_ID) %>%
  filter(WEIGHT_ADMIT <= 300,
         WEIGHT_ADMIT >= 25) %>%
  {n_distinct(.$stay_id)}

## [1] 14923

print("number of valid hourly KDIGO-UO for those stays")

## [1] "number of valid hourly KDIGO-UO for those stays"

kdigo_uo_stage %>% 
  filter(stay_id %in% aki_epi$STAY_ID) %>%
  filter(WEIGHT_ADMIT <= 300,
         WEIGHT_ADMIT >= 25) %>%
  nrow()

## [1] 1558118

 

Inclusion/Exclusion Flowchart

knitr::include_graphics("flow chart.png")

Table 1 - Patient’s characteristics

table_1$SERVICE <- as.factor(table_1$SERVICE)
table_1$admission_age <- as.numeric(table_1$admission_age)
table_1$weight_admit <- as.numeric(table_1$weight_admit)
table_1$height_first <- as.numeric(table_1$height_first)
table_1$creat_peak_72 <- as.numeric(table_1$creat_peak_72)
table_1$creat_peak_72 <- as.numeric(table_1$creat_peak_72)
table_1$creat_last <- as.numeric(table_1$creat_last)

uo_for_table1 <- uo_rate_including_null_collection_period %>%
  drop_na(TIME_INTERVAL) %>%
  group_by(STAY_ID) %>%
  summarise(
    count = n(),
    volumes = mean(VALUE, na.rm = TRUE),
    collection_times = mean(TIME_INTERVAL, na.rm = TRUE),
    rates = mean(HOURLY_RATE, na.rm = TRUE),
    ml_kg_hr = mean(HOURLY_RATE / WEIGHT_ADMIT, na.rm = TRUE)
  )

t1a <- table_1 %>%
  select(
    admission_age,
    gender,
    weight_admit,
    creat_first,
    creat_peak_72,
    creat_last,
    icu_days,
    rrt_binary,
    hospital_expire_flag
  ) %>%
  tbl_summary(
    type = list(
      c(hospital_expire_flag,
        rrt_binary) ~ "dichotomous",
      c(admission_age,
        weight_admit,
        creat_peak_72,
        creat_first,
        creat_last) ~ "continuous"
    ),
    statistic = c(admission_age,
                  weight_admit,
                  creat_first,
                  creat_peak_72,
                  creat_last) ~ "{mean} ({sd})",
    digits = list(icu_days ~ 1),
    missing = "no",
    label = list(
      admission_age ~ "Age at Hospital Admission, years",
      gender ~ "Gender",
      creat_first ~ "First Creatinine in Hospital, mg/dL",
      creat_peak_72 ~ "Peak Creatinine at first days, mg/dL",
      creat_last ~ "Hospital Discharge Creatinine, mg/dL",
      weight_admit ~ "Weight at ICU Admission, kg",
      icu_days ~ "Time in ICU, days",
      rrt_binary ~ "Renal replacement therapy",
      hospital_expire_flag ~ "Hospital Mortality"
    )
  ) %>%
  add_n() 

display_prec <- function(x)
  mean(x) * 100

t1b <- uo_for_table1 %>%
  select(count,
         volumes,
         collection_times,
         rates,
         ml_kg_hr) %>%
  tbl_summary(
    type = list(
      c(volumes,
        collection_times,
        rates,
        ml_kg_hr) ~ "continuous"
    ),
    statistic = list(
      c(volumes,
        collection_times,
        rates,
        ml_kg_hr) ~ "{mean} ({sd})"
    ),
    missing = "no",
    label = list(
      count ~ "Number of Measurements",
      volumes ~ "Average Volumes, mL",
      collection_times ~ "Average Collection Times, minutes",
      rates ~ "Average Rates, mL/hr",
      ml_kg_hr ~ "Average Rate to Weight, mL/hr/kg"
    )
  ) %>%
  add_n() 

tbl_stack(
  list(t1a, t1b),
  group_header = c("ICU Stay", "UO Charting Across ICU Stay")
) %>%
  as_gt() %>%
  tab_source_note(source_note = "The variables age at hospital admission, gender, CCI, CKD, ethnicity, time in hospital, and mortality are measured for each hospital admission and might be counted for more than one ICU stay. All other variables are measured individually for each ICU stay. The variables average collection times, average rates, and average rate to weight are only presented for ICU stays with at least two UO measurements for the same compartment, and the latter also required weight at admission. AKI variables were presented for ICU stays with at least one hour with a non-null KDIGO-UO stage. The variables average AKI duration and total time in AKI were specifically reported for ICU stays with at least one AKI event.") %>%
  tab_source_note(source_note = "AKI: Acute Kidney Injury; CCI: Charlson Comorbidity Index; CKD Stage 1-4: Chronic Kidney Disease excluding end-stage-renal-disease; ICU: Intensive Care Unit; KDIGO: Kidney Disease: Improving Global Outcomes; SOFA: Sequential Organ Failure Assessment; UO: Urine Output.")

Characteristic	N	N = 17,7361
ICU Stay
Age at Hospital Admission, years	17,736	63 (15)
Gender	17,320
F		5,705 (33%)
M		11,615 (67%)
Weight at ICU Admission, kg	17,143	81 (15)
First Creatinine in Hospital, mg/dL	17,677	1.21 (1.06)
Peak Creatinine at first days, mg/dL	16,962	1.34 (1.16)
Hospital Discharge Creatinine, mg/dL	17,677	1.10 (0.84)
Time in ICU, days	17,736	1.3 (0.9, 4.7)
Renal replacement therapy	17,736	1,047 (5.9%)
Hospital Mortality	17,736	2,010 (11%)
UO Charting Across ICU Stay
Number of Measurements	17,736	20 (14, 74)
Average Volumes, mL	17,736	129 (74)
Average Collection Times, minutes	17,736	88 (100)
Average Rates, mL/hr	17,736	106 (97)
Average Rate to Weight, mL/hr/kg	17,143	1.36 (1.26)
The variables age at hospital admission, gender, CCI, CKD, ethnicity, time in hospital, and mortality are measured for each hospital admission and might be counted for more than one ICU stay. All other variables are measured individually for each ICU stay. The variables average collection times, average rates, and average rate to weight are only presented for ICU stays with at least two UO measurements for the same compartment, and the latter also required weight at admission. AKI variables were presented for ICU stays with at least one hour with a non-null KDIGO-UO stage. The variables average AKI duration and total time in AKI were specifically reported for ICU stays with at least one AKI event.
AKI: Acute Kidney Injury; CCI: Charlson Comorbidity Index; CKD Stage 1-4: Chronic Kidney Disease excluding end-stage-renal-disease; ICU: Intensive Care Unit; KDIGO: Kidney Disease: Improving Global Outcomes; SOFA: Sequential Organ Failure Assessment; UO: Urine Output.
1 Mean (SD); n (%); Median (Q1, Q3)

 

Age

Sage_a <- table_1 %>%
  dplyr::summarise(N = n(),
                   Mean = round(mean(admission_age),2),
                   SD = round(sd(admission_age),2),
                   '5th' = round(quantile(admission_age, 0.05),2),
                   '10th' = round(quantile(admission_age, 0.1),2),
                   '25th' = round(quantile(admission_age, 0.25),2),
                   '50th' = round(quantile(admission_age, 0.50),2),
                   '75th' = round(quantile(admission_age, 0.75),2),
                   '95th' = round(quantile(admission_age, 0.95),2),
                   Min = round(min(admission_age),2),
                   Max = round(max(admission_age),2)
  ) %>% gt() %>%
  fmt_number(use_seps = TRUE, decimals = 2)
Sage_a

N	Mean	SD	5th	10th	25th	50th	75th	95th	Min	Max
17,736.00	63.36	15.22	29.00	45.00	55.00	65.00	75.00	85.00	29.00	85.00

Sage_b <- ggplot() + 
  geom_histogram(aes(x = admission_age
                     ), data=table_1, binwidth = 1) + 
  labs(
        x = "Age (years)",
        y = "Frequency"
      )

Sage_b

 

Weight

Sweight_a <- table_1 %>%
  drop_na(weight_admit) %>%
  dplyr::summarise(N = n(),
                   Mean = round(mean(weight_admit),2),
                   SD = round(sd(weight_admit),2),
                   '5th' = round(quantile(weight_admit, 0.05),2),
                   '10th' = round(quantile(weight_admit, 0.1),2),
                   '25th' = round(quantile(weight_admit, 0.25),2),
                   '50th' = round(quantile(weight_admit, 0.50),2),
                   '75th' = round(quantile(weight_admit, 0.75),2),
                   '95th' = round(quantile(weight_admit, 0.95),2),
                   Min = round(min(weight_admit),2),
                   Max = round(max(weight_admit),2)
  ) %>% gt() %>%
  fmt_number(use_seps = TRUE, decimals = 2)
Sweight_a

N	Mean	SD	5th	10th	25th	50th	75th	95th	Min	Max
17,143.00	80.54	14.52	55.00	65.00	75.00	75.00	85.00	105.00	55.00	115.00

Sweight_b <- ggplot() + 
  geom_histogram(aes(x = weight_admit
                     ), data=table_1, binwidth = 5) + 
  labs(
        # title = "Hourly-Adjusted UO per Kilogram",
        x = "Weight (kg)",
        y = "Frequency"
      ) +
  xlim(0, 300)

Sweight_b

 

Table 2 - UO records characteristics

uo_rate$SOURCE <- as.factor(uo_rate$SOURCE)
uo_rate$SERVICE <- as.factor(uo_rate$SERVICE)

uo_rate %>%
  select(VALUE, TIME_INTERVAL, SOURCE, SERVICE) %>%
  tbl_summary(by=SERVICE) %>%
  add_overall()

Characteristic	Overall N = 1,489,5271	IC N = 1,131,2791	IC&MC N = 349,0141	MC&IC N = 9,2341
VALUE	100 (55, 170)	100 (50, 160)	110 (70, 180)	75 (40, 130)
TIME_INTERVAL	60 (60, 120)	60 (60, 120)	60 (60, 120)	60 (60, 120)
SOURCE
Foley	1,458,673 (98%)	1,107,795 (98%)	341,966 (98%)	8,912 (97%)
Ileoconduit	4,749 (0.3%)	2,441 (0.2%)	2,248 (0.6%)	60 (0.6%)
L Nephrostomy	2,404 (0.2%)	2,272 (0.2%)	129 (<0.1%)	3 (<0.1%)
R Nephrostomy	3,423 (0.2%)	2,744 (0.2%)	498 (0.1%)	181 (2.0%)
Suprapubic	11,994 (0.8%)	10,503 (0.9%)	1,446 (0.4%)	45 (0.5%)
Void	8,284 (0.6%)	5,524 (0.5%)	2,727 (0.8%)	33 (0.4%)
1 Median (Q1, Q3); n (%)

---

Raw data analysis

Collection Periods

bS6_a <- uo_rate %>% group_by(SOURCE) %>%
   dplyr::summarise(N = n(),
                   Mean = round(mean(TIME_INTERVAL),0),
                   SD = round(sd(TIME_INTERVAL),0),
                   '5th' = round(quantile(TIME_INTERVAL, 0.05),0),
                   '10th' = round(quantile(TIME_INTERVAL, 0.1),0),
                   '25th' = round(quantile(TIME_INTERVAL, 0.25),0),
                   '50th' = round(quantile(TIME_INTERVAL, 0.50),0),
                   '75th' = round(quantile(TIME_INTERVAL, 0.75),0),
                   '95th' = round(quantile(TIME_INTERVAL, 0.95),0),
                   Min = round(min(TIME_INTERVAL),0),
                   Max = round(max(TIME_INTERVAL),0)
  ) %>% 
  arrange(desc(N)) %>% 
  gt() %>%
  fmt_number(use_seps = TRUE, decimals = 0)

bS6_a

SOURCE	N	Mean	SD	5th	10th	25th	50th	75th	95th	Min	Max
Foley	1,458,673	86	87	60	60	60	60	120	120	1	26,940
Suprapubic	11,994	89	49	60	60	60	60	120	180	1	1,080
Void	8,284	158	383	60	60	60	120	180	420	1	21,420
Ileoconduit	4,749	129	794	60	60	60	120	120	228	1	35,940
R Nephrostomy	3,423	127	326	60	60	60	120	120	300	1	9,720
L Nephrostomy	2,404	139	361	60	60	60	120	120	300	1	9,720

bS6_b <- ggplot(data = uo_rate, aes(x = TIME_INTERVAL / 60)) +
  geom_histogram(binwidth = 1) +
  facet_wrap(~factor(SOURCE, levels=c('Foley', 'Suprapubic', 'Ileoconduit',
                         'Void',
                         'R Nephrostomy', 'L Nephrostomy')), scales = "free") +
  xlim(-1, 20) +
  labs(
          x = "Time interval (hr)",
          y = "Frequency"
        )

bS6_b

 

Volumes and Collection Periods

bS7_a <- uo_rate %>% group_by(SOURCE) %>%
  dplyr::summarise(N = n(),
                   Mean = round(mean(VALUE),0),
                   SD = round(sd(VALUE),0),
                   '5th' = round(quantile(VALUE, 0.05),0),
                   '10th' = round(quantile(VALUE, 0.1),0),
                   '25th' = round(quantile(VALUE, 0.25),0),
                   '50th' = round(quantile(VALUE, 0.50),0),
                   '75th' = round(quantile(VALUE, 0.75),0),
                   '95th' = round(quantile(VALUE, 0.95),0),
                   Min = round(min(VALUE),0),
                   Max = round(max(VALUE),0)
  ) %>% 
  arrange(desc(N)) %>% 
  gt() %>%
  fmt_number(use_seps = TRUE, decimals = 0)

bS7_a

SOURCE	N	Mean	SD	5th	10th	25th	50th	75th	95th	Min	Max
Foley	1,458,673	131	113	20	30	55	100	170	360	0	4,450
Suprapubic	11,994	124	118	10	20	50	90	160	370	0	2,000
Void	8,284	213	179	20	45	90	160	300	550	0	1,600
Ileoconduit	4,749	130	114	20	30	60	100	160	350	0	1,400
R Nephrostomy	3,423	65	85	0	0	15	40	85	200	0	1,200
L Nephrostomy	2,404	73	92	0	1	15	40	100	260	0	900

bS7_b <- ggplot(data = uo_rate, aes(x = VALUE)) +
  facet_wrap(~factor(SOURCE, levels=c('Foley', 'Suprapubic', 'Ileoconduit',
                         'Void',
                         'R Nephrostomy', 'L Nephrostomy')), scales = "free") +  geom_histogram(binwidth = 50) +
  xlim(-25, 1100) +
  labs(
        # title = "Volumes",
        x = "Volume (ml)",
        y = "Frequency"
      ) +
      theme(
        plot.title = element_text(color = "#0099F8", size = 16, face = "bold"),
        plot.subtitle = element_text(size = 10, face = "bold"),
        plot.caption = element_text(face = "italic")
      )

bS7_b

 

Records of zero volume

The proportion of zero value UO measurements:

uo_rate_count <- uo_rate %>% 
  count(SOURCE, sort = TRUE)

uo_rate_0_count <- uo_rate %>% 
  filter(VALUE == 0) %>% 
  count(SOURCE, sort = TRUE)

count_uo_zero_vs_all <- left_join(uo_rate_count, 
                                  uo_rate_0_count, by = "SOURCE") %>% 
  mutate(PROPORTION = n.y / n.x) %>%
  pivot_longer(cols = n.y:n.x, names_to = "type")
  
bS7_c <- count_uo_zero_vs_all %>%
  ggplot(aes(x=reorder(SOURCE, -value), y=value, fill=type)) +
    geom_bar(position="fill", stat="identity") +
    xlab("") +
    ylab("") +
    scale_fill_brewer(palette="Paired") +  
    geom_text(aes(label=ifelse(type == "n.y", paste0((round(PROPORTION, 3) * 100), "%"), "")), 
          color="black", 
          size=3, 
          vjust=-1,
          position="fill") +
    theme_minimal() +
    theme(axis.text.y=element_blank()) +
    theme(legend.position="none") +
  theme(axis.text.y=element_blank()) +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))+
     labs(
         # title = "Proportion of zero value raw output count",
       )
bS7_c

uo_rate_0 <- uo_rate %>% filter(VALUE == 0) 
bS7_d <- uo_rate_0 %>% group_by(SOURCE) %>%
  dplyr::summarise(N = n(),
                   Mean = round(mean(TIME_INTERVAL),0),
                   SD = round(sd(TIME_INTERVAL),0),
                   '5th' = round(quantile(TIME_INTERVAL, 0.05),0),
                   '10th' = round(quantile(TIME_INTERVAL, 0.1),0),
                   '25th' = round(quantile(TIME_INTERVAL, 0.25),0),
                   '50th' = round(quantile(TIME_INTERVAL, 0.50),0),
                   '75th' = round(quantile(TIME_INTERVAL, 0.75),0),
                   '95th' = round(quantile(TIME_INTERVAL, 0.95),0),
                   Min = round(min(TIME_INTERVAL),0),
                   Max = round(max(TIME_INTERVAL),0)
  ) %>% 
  arrange(desc(N)) %>% 
  gt() %>%
  fmt_number(use_seps = TRUE, decimals = 0)

bS7_d

SOURCE	N	Mean	SD	5th	10th	25th	50th	75th	95th	Min	Max
Foley	22,280	114	245	60	60	60	60	120	300	1	20,868
R Nephrostomy	355	163	527	60	60	60	120	120	360	1	9,720
Suprapubic	337	106	87	60	60	60	60	120	240	1	600
Void	265	257	475	60	60	60	120	240	996	1	4,560
L Nephrostomy	231	184	373	60	60	60	120	120	540	59	5,040
Ileoconduit	65	587	3,425	60	60	60	60	120	456	1	27,480

---

Adjusting for hourly UO

UO Rate

bS8_a <- uo_rate %>% group_by(SOURCE) %>%
    dplyr::summarise(N = n(),
                   Mean = round(mean(HOURLY_RATE),0),
                   SD = round(sd(HOURLY_RATE),0),
                   '5th' = round(quantile(HOURLY_RATE, 0.05),0),
                   '10th' = round(quantile(HOURLY_RATE, 0.1),0),
                   '25th' = round(quantile(HOURLY_RATE, 0.25),0),
                   '50th' = round(quantile(HOURLY_RATE, 0.50),0),
                   '75th' = round(quantile(HOURLY_RATE, 0.75),0),
                   '95th' = round(quantile(HOURLY_RATE, 0.95),0),
                   Min = round(min(HOURLY_RATE),0),
                   Max = round(max(HOURLY_RATE),0)
  ) %>% 
  arrange(desc(N)) %>% 
  gt() %>%
  fmt_number(use_seps = TRUE, decimals = 0)

bS8_a

SOURCE	N	Mean	SD	5th	10th	25th	50th	75th	95th	Min	Max
Foley	1,458,673	105	235	13	25	50	75	125	270	0	42,000
Suprapubic	11,994	100	389	5	17	38	65	120	265	0	30,000
Void	8,284	119	264	10	25	50	90	150	300	0	21,000
Ileoconduit	4,749	83	125	10	25	40	65	100	208	0	7,200
R Nephrostomy	3,423	43	82	0	0	9	30	60	125	0	3,600
L Nephrostomy	2,404	52	103	0	1	9	30	60	200	0	3,000

bS8_b <- ggplot(data = uo_rate, aes(x = HOURLY_RATE)) +
  geom_histogram(binwidth = 20) +
  facet_wrap(~factor(SOURCE, levels=c('Foley', 'Suprapubic', 'Ileoconduit',
                         'Void',
                         'R Nephrostomy', 'L Nephrostomy')), scales = "free") +  xlim(-10, 500) +
  labs(
        # title = "UO Rates",
        # subtitle = "by source",
        x = "Rate (ml/hr)",
        y = "Frequency"
      ) +
      theme(
        plot.title = element_text(color = "#0099F8", size = 16, face = "bold"),
        plot.subtitle = element_text(size = 10, face = "bold"),
        plot.caption = element_text(face = "italic")
      )

bS8_b

 

Low UO Rate Analysis

The association between UO rates and collection periods, smoothed conditional means for records of Foley catheter.

bS8_d <- uo_rate %>%
  filter(SOURCE == "Foley",
         HOURLY_RATE < 500) %>%
  # slice_head(n = 50000) %>%
  ggplot(aes(x = HOURLY_RATE, y = TIME_INTERVAL)) +
  geom_smooth(se = TRUE, alpha = 0.1, linewidth = 1) +
  # geom_smooth(
  #   se = FALSE,
  #   method = "lm",
  #   linetype = "dashed",
  #   color = "red", 
  #   linewidth = 0.3
  # ) +
  geom_hline(yintercept = 60,
             size = 0.3,
             color = "#cccccc") +
  geom_vline(
    xintercept = 20,
    size = 0.3,
    color = "black",
    linetype = "dotdash"
  ) +
  scale_x_continuous(breaks = c(0, 20, 50, 100, 200, 300, 400, 500)) +
  scale_y_continuous(breaks = c(0, 60, 100, 200)) +
  coord_cartesian(xlim = c(0, 500), ylim = c(0, 200)) +
  labs(x = "Urine output rate (ml/hr)", y = "Collection periods (min)") +
  theme_classic()

bS8_d

Quantile analysis for collection periods as a function of rates.

uo_rate_qreg <- uo_rate %>%
  filter(SOURCE == "Foley",
         HOURLY_RATE < 500) %>%
  arrange(desc(CHARTTIME)) %>%
  slice_head(n = 500000)

# #### Quantile
quantile_reg <- rq(TIME_INTERVAL ~
                     HOURLY_RATE,
                   seq(0.10, 0.90, by = 0.10),
                   # c(.05, .1, .25, .5, .75, .90, .95),
                   data = uo_rate_qreg)

# summary(quantile_reg, se = "iid") %>% 
#   plot()

### OLS
lm <- lm(data=uo_rate_qreg,
         formula =  TIME_INTERVAL ~
           HOURLY_RATE)

ols <- as.data.frame(coef(lm))
ols.ci <- as.data.frame(confint(lm, level = 0.95))
ols2 <- cbind(ols, ols.ci)
ols2 <- tibble::rownames_to_column(ols2, var="term")



#### Quantile
bS8_e <- quantile_reg %>%
  tidy(se.type = "iid", conf.int = TRUE, conf.level = 0.95) %>%
  filter(!grepl("factor", term)) %>%
  ggplot(aes(x=tau,y=estimate)) +
  theme_classic() +
  theme(
    strip.background = element_blank(),
    #strip.text.x = element_blank()
  ) +
  scale_y_continuous(limits = symmetric_limits) +
  scale_x_continuous(breaks = scales::pretty_breaks(n = 12)) +
  ##### quantilie results
  geom_point(color="#27408b", size = 0.3)+ 
  geom_line(color="black", linetype = "dotdash", size = 0.3)+ 
  geom_ribbon(aes(ymin=conf.low,ymax=conf.high),alpha=0.25, fill="#555555")+
  facet_wrap(~term, scales="free", ncol=1)+
  ##### OLS results
  geom_hline(data = ols2, aes(yintercept= `coef(lm)`), lty=1, color="red", size=0.3)+
  geom_hline(data = ols2, aes(yintercept= `2.5 %`), lty=2, color="red", size=0.3)+
  geom_hline(data = ols2, aes(yintercept= `97.5 %`), lty=2, color="red", size=0.3)+
  #### Lines
   geom_hline(yintercept = 0, size=0.3) 

bS8_e

# Visualization for Quantile Regression with some tau values: 
intercept_slope <- quantile_reg %>% 
  coef() %>% 
  t() %>% 
  data.frame() %>% 
  rename(intercept = X.Intercept., slope = HOURLY_RATE) %>% 
  mutate(quantile = row.names(.))


bS8_f <-
  ggplot() +
  geom_point(data = uo_rate_qreg, aes(HOURLY_RATE, TIME_INTERVAL),
    alpha = 0.2,
    size = 0.5,
    stroke = 0.5,
    width = 2,
    height = 2
  ) +
  geom_abline(data = intercept_slope, aes(
    intercept = intercept,
    slope = slope,
    color = quantile
  ),
  linewidth=1) +
  theme_minimal() +
  labs(x = "Urine output rate (ml/hr)", y = "Collection periods (min)") +
  coord_cartesian(xlim = c(0, 500), ylim = c(0, 500)) +
  theme(legend.position="bottom")

bS8_f

 

Collection periods for UO rate 20ml/hr or below

uo_rate %>% 
  filter(HOURLY_RATE <= 20) %>%
  group_by(SOURCE) %>%
   dplyr::summarise(N = n(),
                   Mean = round(mean(TIME_INTERVAL),0),
                   SD = round(sd(TIME_INTERVAL),0),
                   '5th' = round(quantile(TIME_INTERVAL, 0.05),0),
                   '10th' = round(quantile(TIME_INTERVAL, 0.1),0),
                   '25th' = round(quantile(TIME_INTERVAL, 0.25),0),
                   '50th' = round(quantile(TIME_INTERVAL, 0.50),0),
                   '75th' = round(quantile(TIME_INTERVAL, 0.75),0),
                   '95th' = round(quantile(TIME_INTERVAL, 0.95),0),
                   Min = round(min(TIME_INTERVAL),0),
                   Max = round(max(TIME_INTERVAL),0)
  ) %>% 
  arrange(desc(N)) %>% 
  gt() %>%
  fmt_number(use_seps = TRUE, decimals = 0)

SOURCE	N	Mean	SD	5th	10th	25th	50th	75th	95th	Min	Max
Foley	123,041	123	270	60	60	60	60	120	300	1	26,940
Suprapubic	1,588	113	88	60	60	60	116	120	240	1	1,080
R Nephrostomy	1,500	168	485	60	60	60	120	120	360	1	9,720
L Nephrostomy	1,073	194	532	60	60	60	120	180	480	59	9,720
Void	697	386	1,243	60	60	60	120	360	1,032	1	21,420
Ileoconduit	450	397	2,561	60	60	60	120	120	480	1	35,940

uo_rate %>% 
  filter(HOURLY_RATE <= 20) %>%
ggplot(aes(x = TIME_INTERVAL / 60)) +
  geom_histogram(binwidth = 1) +
  facet_wrap(~factor(SOURCE, levels=c('Foley', 'Suprapubic', 'Ileoconduit',
                         'Void',
                         'R Nephrostomy', 'L Nephrostomy')), scales = "free") +  xlim(-1, 20) +
  labs(
          title = "Collection periods for UO rate 20ml/hr or below",
          subtitle = "by source",
          x = "Time interval (hr)",
          y = "Frequency"
        ) +
        theme(
          plot.title = element_text(color = "#0099F8", size = 16, face = "bold"),
          plot.subtitle = element_text(size = 10, face = "bold"),
          plot.caption = element_text(face = "italic")
        )

 

Mean Rate

Mean UO rate weighted by tyme and grouped by source:

bS8_c <- uo_rate %>% 
  group_by(SOURCE) %>%
  summarise(weighted_mean_rate = weighted.mean(HOURLY_RATE, TIME_INTERVAL)) %>%
  gt() %>%
  fmt_number(use_seps = TRUE, decimals = 1) %>%
  cols_label(
    SOURCE = "Source",
    weighted_mean_rate = "Weighted mean rate (ml/hr)"
  )

bS8_c

Source	Weighted mean rate (ml/hr)
Foley	91.2
Ileoconduit	60.5
L Nephrostomy	31.8
R Nephrostomy	30.7
Suprapubic	83.3
Void	81.1

 

Hourly-adjusted UO

# 2065999
bS9_a <- hourly_uo %>% drop_na(HOURLY_WEIGHTED_MEAN_RATE) %>%
    dplyr::summarise(N = n(),
                   Mean = round(mean(HOURLY_WEIGHTED_MEAN_RATE),0),
                   SD = round(sd(HOURLY_WEIGHTED_MEAN_RATE),0),
                   '5th' = round(quantile(HOURLY_WEIGHTED_MEAN_RATE, 0.05),0),
                   '10th' = round(quantile(HOURLY_WEIGHTED_MEAN_RATE, 0.1),0),
                   '25th' = round(quantile(HOURLY_WEIGHTED_MEAN_RATE, 0.25),0),
                   '50th' = round(quantile(HOURLY_WEIGHTED_MEAN_RATE, 0.50),0),
                   '75th' = round(quantile(HOURLY_WEIGHTED_MEAN_RATE, 0.75),0),
                   '95th' = round(quantile(HOURLY_WEIGHTED_MEAN_RATE, 0.95),0),
                   Min = round(min(HOURLY_WEIGHTED_MEAN_RATE),0),
                   Max = round(max(HOURLY_WEIGHTED_MEAN_RATE),0)
  ) %>% 
  arrange(desc(N)) %>% 
  gt() %>%
  fmt_number(use_seps = TRUE, decimals = 0)

bS9_a

N	Mean	SD	5th	10th	25th	50th	75th	95th	Min	Max
2,140,972	91	79	5	19	43	72	117	238	0	3,720

bS9_b <- hourly_uo %>% drop_na(HOURLY_WEIGHTED_MEAN_RATE) %>%
ggplot(aes(x = HOURLY_WEIGHTED_MEAN_RATE)) +
  geom_histogram(binwidth = 25) +
  xlim(-10, 500) + 
  labs(
        # title = "Hourly-Adjusted UO",
        x = "Hourly UO (ml)",
        y = "Frequency"
      )

bS9_b

 

Simple Sum Comparison

Showing proportion of hours with less than 100ml difference):

adj_uo_diff <- hourly_uo %>%
  select(HOURLY_WEIGHTED_MEAN_RATE, SIMPLE_SUM) %>%
  filter(!is.na(HOURLY_WEIGHTED_MEAN_RATE)) %>%
  mutate(hourly_diff = abs(HOURLY_WEIGHTED_MEAN_RATE - SIMPLE_SUM)) %>%
  mutate(
    cutoff_10 = if_else(hourly_diff < 10, 1, 0),
    cutoff_50 = if_else(hourly_diff < 50, 1, 0),
    cutoff_100 = if_else(hourly_diff < 100, 1, 0),
    cutoff_150 = if_else(hourly_diff < 150, 1, 0),
    cutoff_200 = if_else(hourly_diff < 200, 1, 0)
  )

my_order <- c("<10", "<50", "<100", "<150", "<200")

bS11 <- adj_uo_diff %>%
  select(cutoff_10,
         cutoff_50,
         cutoff_100,
         cutoff_150,
         cutoff_200) %>%
  pivot_longer(cols = contains("cutoff")) %>%
  transmute(name = case_when(
    name == "cutoff_10" ~ "<10",
    name == "cutoff_50" ~ "<50",
    name == "cutoff_100" ~ "<100",
    name == "cutoff_150" ~ "<150",
    name == "cutoff_200" ~ "<200"
  ),
  value) %>%
  group_by(name) %>%
  summarise(agreement = paste0(round(mean(value) * 100, 1), "%"),
            non_agreement = paste0(round((1 - mean(
              value
            )) * 100, 1), "%")) %>%
  arrange(match(name, my_order)) %>%
  gt() %>%
  # tab_header(
  #   title = md("**Comparison of Hourly-Adjusted UO and Simple Summation**"),
  # ) %>%
  cols_label(
    name = "Cut-off (ml)",
    agreement = "Proportion  of Agreement",
    non_agreement = "Proportion  of Disagreement"
  ) %>%
  cols_align(
    align = "center"
  ) %>%
  tab_source_note(source_note = "The table demonstrates the significance of hourly adjustment for accuracy by presenting the variance between the adjusted values and the simple hourly summation. Cut-off values are based on the absolute difference between the hourly-adjusted UO and a simple hourly summation of UO. Measurements charted on the hour were included with the previous time interval. ")

adj_uo_diff <- hourly_uo %>%
  select(HOURLY_WEIGHTED_MEAN_RATE, SIMPLE_SUM) %>%
  filter(!is.na(HOURLY_WEIGHTED_MEAN_RATE)) %>%
  mutate(no_diff = 
           ifelse((is.na(HOURLY_WEIGHTED_MEAN_RATE) &
                  is.na(SIMPLE_SUM)) |
             (!is.na(HOURLY_WEIGHTED_MEAN_RATE) &
                  !is.na(SIMPLE_SUM) &
                    abs(HOURLY_WEIGHTED_MEAN_RATE-SIMPLE_SUM) < 100), 
                  1, 
                  0),
         .keep = "none")

bS11

Cut-off (ml)	Proportion of Agreement	Proportion of Disagreement
<10	29.5%	70.5%
<50	60.1%	39.9%
<100	82.2%	17.8%
<150	92%	8%
<200	96.1%	3.9%
The table demonstrates the significance of hourly adjustment for accuracy by presenting the variance between the adjusted values and the simple hourly summation. Cut-off values are based on the absolute difference between the hourly-adjusted UO and a simple hourly summation of UO. Measurements charted on the hour were included with the previous time interval.

mean(adj_uo_diff$no_diff)

## [1] 0.8224012

 

Hourly UO Per Kilogram

# 1993167
bS9_c <- uo_ml_kg_hr %>%
  filter(WEIGHT_ADMIT <= 300,
         WEIGHT_ADMIT >= 25) %>%
  dplyr::summarise(N = n(),
                   Mean = round(mean(ML_KG_HR),2),
                   SD = round(sd(ML_KG_HR),2),
                   '5th' = round(quantile(ML_KG_HR, 0.05),2),
                   '10th' = round(quantile(ML_KG_HR, 0.1),2),
                   '25th' = round(quantile(ML_KG_HR, 0.25),2),
                   '50th' = round(quantile(ML_KG_HR, 0.50),2),
                   '75th' = round(quantile(ML_KG_HR, 0.75),2),
                   '95th' = round(quantile(ML_KG_HR, 0.95),2),
                   Min = round(min(ML_KG_HR),2),
                   Max = round(max(ML_KG_HR),2)
  ) %>% 
  gt() %>%
  fmt_number(use_seps = TRUE, decimals = 2)

bS9_c

N	Mean	SD	5th	10th	25th	50th	75th	95th	Min	Max
2,063,720.00	1.18	1.06	0.07	0.24	0.54	0.92	1.50	3.09	0.00	67.64

mean_log <- log(mean(uo_ml_kg_hr$ML_KG_HR))
sd_log <- log(sd(uo_ml_kg_hr$ML_KG_HR))
bS9_d <- ggplot() + 
  xlim(0, 2) + 
  geom_histogram(aes(x = ML_KG_HR
                     # , y =..density..
                     ), data=(uo_ml_kg_hr %>%
                                filter(WEIGHT_ADMIT <= 300,
                                       WEIGHT_ADMIT >= 25)), 
                 binwidth = 0.02) + 
  # stat_function(fun = dlnorm, args = list(meanlog = mean_log, sdlog = sd_log, log = FALSE), size=1, color='gray') +
  labs(
        # title = "Hourly-Adjusted UO per Kilogram",
        x = "Hourly volume to kg (ml/hr/kg)",
        y = "Frequency"
      )

bS9_d

KDIGO Criteria Avarage-UO and Consecutive-UO Compariton

aumc_kdigo_inter_aki_table <- akis_all_long %>%
  filter(group != "old") %>%
  drop_na(prevalnce_admit) %>%
  transmute(
    group = case_when(
      group == "newcons" ~ 'UO-Consecutive',
      group == "newmean" ~ 'UO-Average',
      group == "old" ~ 'Block summation'
    ),
    aki_binary = if_else(max_stage > 0, 1, 0),
    max_stage = if_else(max_stage == 0, NA, max_stage),
    prevalnce_admit
  ) %>%
  tbl_summary(
    by = "group",
    missing = "no",
    digits = everything() ~ c(0, 1),
    label = list(
      aki_binary ~ "Oliguric-AKI on the first days",
      prevalnce_admit ~ "Prevalence at admission",
      max_stage ~ "Maximum KDIGO staging"
    )
  )  %>%
  modify_column_indent(columns = label, rows = c(FALSE, TRUE)) %>%
  modify_column_indent(
    columns = label,
    rows = c(FALSE, FALSE, TRUE, TRUE, TRUE),
    double_indent = TRUE
  ) %>%
  add_p()

aumc_kdigo_inter_aki_table

Characteristic	UO-Average N = 14,9231	UO-Consecutive N = 14,9231	p-value2
Oliguric-AKI on the first days	7,429 (49.8%)	4,688 (31.4%)	<0.001
Maximum KDIGO staging			<0.001
1	2,509 (33.8%)	2,456 (52.4%)
2	3,741 (50.4%)	1,589 (33.9%)
3	1,179 (15.9%)	643 (13.7%)
Prevalence at admission	1,962 (13.1%)	1,024 (6.86%)	<0.001
1 n (%)
2 Pearson’s Chi-squared test

check for mortality overdispration:

model_binom <- glm(
  mortality_30 ~ stage_newcons,
  family = binomial,
  data = (
    akis_all_long %>%
      filter(group == "newcons") %>%
      mutate(stage_newcons = as.factor(max_stage))
  )
)

model_overdispersed <- glm(
  mortality_30 ~ stage_newcons,
  family = quasibinomial,
  data = (
    akis_all_long %>%
      filter(group == "newcons") %>%
      mutate(stage_newcons = as.factor(max_stage))
  )
)

pchisq(summary(model_overdispersed)$dispersion * model_binom$df.residual, 
       model_binom$df.residual, lower = F)

## [1] 0.4886805

aumc_kdigo_inter_cons_mean <- akis_all_long %>%
  transmute(STAY_ID,
            group,
            first_stage,
            max_stage,
            aki_above_2 = max_stage > 1) %>%
  filter(aki_above_2 == TRUE,
         (group != "old")) %>%
  mutate(group = case_when(
    group == "newcons" ~ 'UO-Consecutive',
    group == "newmean" ~ 'UO-Average',
    group == "old" ~ 'Block summation'
  )) %>%
  left_join(table_1, by = "STAY_ID") %>%
  select(
    group,
    admission_age,
    weight_admit,
    gender,
    creat_first,
    creat_peak_72,
    creat_last,
  # scr_baseline,
    icu_days,
    rrt_binary,
    hospital_expire_flag
  ) %>%
  tbl_summary(
    by = group,
    type = list(
      c(rrt_binary) ~ "dichotomous",
      c(admission_age,
        weight_admit,
        creat_first,
        creat_peak_72,
        creat_last) ~ "continuous"
    ),
    statistic = c(admission_age,
                  weight_admit,
                  creat_first,
                  creat_peak_72,
                  creat_last) ~ "{mean} ({sd})",
    missing = "no",
    label = list(
      admission_age ~ "Age at Hospital Admission, years",
      gender ~ "Gender",
      weight_admit ~ "Weight at ICU Admission, kg",
      creat_first ~ "First Creatinine in Hospital, mg/dL",
      creat_peak_72 ~ "Peak Creatinine at first days, mg/dL",
      creat_last ~ "Hospital Discharge Creatinine, mg/dL",
      icu_days ~ "Time in ICU, days",
      rrt_binary ~ "Renal replacement therapy",
      hospital_expire_flag ~ "Hospital Mortality"
    )
  ) %>%
  add_p() %>%
  add_stat(
    fns = everything() ~ add_by_n
  ) %>%
  modify_header(starts_with("add_n_stat") ~ "**N**") %>%
  modify_table_body(
    ~ .x %>%
      dplyr::relocate(add_n_stat_1, .before = stat_1) %>%
      dplyr::relocate(add_n_stat_2, .before = stat_2)
  )

aumc_kdigo_inter_cons_mean

Characteristic	N	UO-Average N = 4,9201	N	UO-Consecutive N = 2,2321	p-value2
Age at Hospital Admission, years	4,920	66 (15)	2,232	67 (14)	0.002
Weight at ICU Admission, kg	4,920	84 (15)	2,232	84 (15)	0.057
Gender	4,845		2,203		0.8
F		1,718 (35%)		775 (35%)
M		3,127 (65%)		1,428 (65%)
First Creatinine in Hospital, mg/dL	4,914	1.40 (1.38)	2,230	1.72 (1.86)	<0.001
Peak Creatinine at first days, mg/dL	4,882	1.79 (1.53)	2,218	2.45 (1.90)	<0.001
Hospital Discharge Creatinine, mg/dL	4,914	1.31 (1.08)	2,230	1.67 (1.35)	<0.001
Time in ICU, days	4,920	4 (2, 10)	2,232	5 (2, 13)	<0.001
Renal replacement therapy	4,920	745 (15%)	2,232	696 (31%)	<0.001
Hospital Mortality	4,920	963 (20%)	2,232	697 (31%)	<0.001
1 Mean (SD); n (%); Median (Q1, Q3)
2 Wilcoxon rank sum test; Pearson’s Chi-squared test

akis_all_long %>%
  group_by(group, max_stage) %>%
  summarise(Propo = sum(mortality_7, na.rm = T) / n()) %>%
  ggplot(aes(max_stage, Propo, fill = group)) + geom_col(position = 'dodge')

aumc_cdplod_mean <- cdplot(as.factor(mortality_7)~max_stage, (akis_all_long %>% 
                                        filter(group == "newmean")),
       col=c("lightgoldenrod", "lightcyan"), 
       ylab = "7 days mortality", xlab ="Max KDIGO-UO stage", main = "CD Plot for UOmean")

aumc_cdplod_mean

## $`1`
## function (v) 
## .approxfun(x, y, v, method, yleft, yright, f, na.rm)
## <bytecode: 0x133ecbf30>
## <environment: 0x133e6ab80>

aumc_cdplod_cons <- cdplot(as.factor(mortality_7)~max_stage, (akis_all_long %>% 
                                        filter(group == "newcons")),
       col=c("lightgoldenrod", "lightcyan"), 
       ylab = "7 days mortality", xlab ="Max KDIGO-UO stage", main = "CD Plot for UOcons")

aumc_cdplod_cons

## $`1`
## function (v) 
## .approxfun(x, y, v, method, yleft, yright, f, na.rm)
## <bytecode: 0x133ecbf30>
## <environment: 0x1143f9d98>

model_mean <- glm(
  mortality_7 ~ MAX_STAGE_NEW_MEAN + FIRST_STAGE_NEW_MEAN,
  data = akis_all_wide,
  family = binomial
)

model_cons <- glm(
  mortality_7 ~ MAX_STAGE_NEW_CONS + FIRST_STAGE_NEW_CONS,
  data = akis_all_wide,
  family = binomial
)

# anova(model_old, model_new, test = 'Chisq')
aumc_kdigo_inter_bic <- BIC(model_mean, model_cons)

aumc_kdigo_inter_bic <- cbind(Model = rownames(aumc_kdigo_inter_bic), aumc_kdigo_inter_bic) %>%
  gt()

aumc_kdigo_inter_bic

Model	df	BIC
model_mean	3	6926.717
model_cons	3	6903.629

akis_all_long_complete <- akis_all_long %>% drop_na(max_stage)

# descriptive
descriptive_tbl <- akis_all_long_complete %>%
  group_by(group, max_stage) %>%
  summarise(
    n = n(),
    dead = sum(mortality_30),
    mortality_prop = sum(mortality_30) / n()
  ) %>%
  group_by(group) %>%
  transmute(
    max_stage,
    "Patients, No. (%)" = paste0(n, " (", round((n / sum(
      n
    )), 2), ")"),
    "Mortality, No. (%)" = paste0(dead, " (", round(mortality_prop, 2), ")")
  ) %>%
  gt()

# glm
m1 <- glm(
  mortality_30 ~ stage_newcons,
  family = binomial,
  data = (
    akis_all_long_complete %>%
      filter(group == "newcons") %>%
      mutate(stage_newcons = as.factor(max_stage))
  )
)

m2 <- glm(
  mortality_30 ~ stage_newmean,
  family = binomial,
  data = (
    akis_all_long_complete %>%
      filter(group == "newmean") %>%
      mutate(stage_newmean = as.factor(max_stage))
  )
)

glm_tbl <- tbl_stack(list(
  tbl_regression(m1, exponentiate = TRUE),
  tbl_regression(m2, exponentiate = TRUE)
)) %>%
  as_gt()

# join tables
glm_tbl_data <- glm_tbl$`_data` %>%
  filter(!is.na(term)) %>%
  transmute(
    group = case_when(
      variable == "stage_newmean" ~ "newmean",
      variable == "stage_newcons" ~ "newcons"
    ),
    max_stage = label,
    "OR (95% CI)" = if_else(label == "0", "1 [Reference]", paste0(round(estimate, 2), " (", ci, ")"), ),
    "P value" = case_when(
      is.na(p.value) ~ "NA",
      p.value < 0.001 ~ "<.001",
      .default = as.character(round(p.value, 2))
    )
  )

descriptive_tbl_data <- descriptive_tbl$`_data` %>%
  mutate(max_stage = as.character(max_stage))

rr_table_data <- left_join(descriptive_tbl_data, glm_tbl_data, by = c("group", "max_stage"))

aumc_kdigo_inter_survival_table <- rr_table_data %>%
  mutate(
    group = case_when(
      group == "newcons" ~ 'UO-Consecutive',
      group == "newmean" ~ 'UO-Average'
    )
  ) %>%
  gt(
    rowname_col = "max_stage",
    groupname_col = "group",
    row_group_as_column = TRUE
  ) %>%
  tab_stubhead(label = "Criteria / Stage") %>%
  tab_spanner(label = "Univariate analysis", columns = c("OR (95% CI)", "P value"))

aumc_kdigo_inter_survival_table

Criteria / Stage		Patients, No. (%)	Mortality, No. (%)	Univariate analysis
				OR (95% CI)	P value
UO-Consecutive	0	10235 (0.69)	605 (0.06)	1 [Reference]	NA
	1	2456 (0.16)	384 (0.16)	2.95 (2.57, 3.38)	<.001
	2	1589 (0.11)	428 (0.27)	5.87 (5.11, 6.73)	<.001
	3	643 (0.04)	312 (0.49)	15 (12.6, 17.9)	<.001
UO-Average	0	7494 (0.5)	363 (0.05)	1 [Reference]	NA
	1	2509 (0.17)	278 (0.11)	2.45 (2.08, 2.88)	<.001
	2	3741 (0.25)	570 (0.15)	3.53 (3.08, 4.06)	<.001
	3	1179 (0.08)	518 (0.44)	15.39 (13.2, 18.0)	<.001

adjusted models:

m1_adj1 <- glm(
  mortality_30 ~ stage_newcons * (
    admission_age + gender + weight_admit + first_stage
  ),
  family = binomial,
  data = (
    akis_all_long_complete %>%
      filter(group == "newcons") %>%
      mutate(stage_newcons = as.factor(max_stage))
  )
)

mrg_efct_m1_adj1 <- avg_comparisons(
  m1_adj1,
  variables = "stage_newcons",
  comparison = "lnoravg",
  transform = exp
)

rr_table_data_m1_adj1 <- mrg_efct_m1_adj1 %>%
    transmute(
      group = "newcons",
      max_stage = as.factor(dplyr::row_number()),
      estimate = paste0(round(estimate, 2), " (",round(conf.low, 2), "-", round(conf.high, 2), ")"),
      p.value
    ) %>% add_row(group = "newcons", max_stage = "0", .before = 0) %>%
    transmute(
    group,
    max_stage,
    or.adj1 = if_else(is.na(estimate), "1 [Reference]", estimate),
    p.adj1 = case_when(
      is.na(p.value) ~ "NA",
      p.value < 0.001 ~ "<.001",
      .default = as.character(round(p.value, 2))
    )
  )

m2_adj1 <- glm(
  mortality_30 ~ stage_newmean * (
    admission_age + weight_admit + gender + first_stage
  ),
  family = binomial,
  data = (
    akis_all_long_complete %>%
      filter(group == "newmean") %>%
      mutate(stage_newmean = as.factor(max_stage))
  )
)

mrg_efct_m2_adj1 <- avg_comparisons(
  m2_adj1,
  variables = "stage_newmean",
  comparison = "lnoravg",
  transform = exp
)

rr_table_data_m2_adj1 <- mrg_efct_m2_adj1 %>%
    transmute(
      group = "newmean",
      max_stage = as.factor(dplyr::row_number()),
      estimate = paste0(round(estimate, 2), " (",round(conf.low, 2), "-", round(conf.high, 2), ")"),
      p.value
    ) %>% add_row(group = "newmean", max_stage = "0", .before = 0) %>%
    transmute(
    group,
    max_stage,
    or.adj1 = if_else(is.na(estimate), "1 [Reference]", estimate),
    p.adj1 = case_when(
      is.na(p.value) ~ "NA",
      p.value < 0.001 ~ "<.001",
      .default = as.character(round(p.value, 2))
    )
  )


aumc_kdigo_inter_survival_table_adj <- rr_table_data %>% 
  left_join(bind_rows(rr_table_data_m1_adj1, rr_table_data_m2_adj1)) %>%
  mutate(
    group = case_when(
      group == "newcons" ~ 'UO-Consecutive',
      group == "newmean" ~ 'UO-Average',
    )
  ) %>%
  gt(
    rowname_col = "max_stage",
    groupname_col = "group",
    row_group_as_column = TRUE
  ) %>%
  cols_label(
    or.adj1 = "OR (95% CI)",
    p.adj1 = "P value"
  ) %>%
  tab_stubhead(label = "Criteria / Stage") %>%
  tab_spanner(label = "Unadjusted OR", columns = c("OR (95% CI)", "P value")) %>%
  tab_spanner(label = "Adjusted Model", columns = c(or.adj1, p.adj1), id = "adj1") %>%
  tab_footnote(
    footnote = "Model include age, weight, gender and whether diagnosed on admission",
    locations = cells_column_spanners(spanners = "adj1")
  ) %>% tab_source_note(source_note = md(
    "All covariates in the adjusted model were significant except for diagnosis at admission for block summation model."
  ))

summaries for all tables:

summary(m1_adj1)

## 
## Call:
## glm(formula = mortality_30 ~ stage_newcons * (admission_age + 
##     gender + weight_admit + first_stage), family = binomial, 
##     data = (akis_all_long_complete %>% filter(group == "newcons") %>% 
##         mutate(stage_newcons = as.factor(max_stage))))
## 
## Coefficients: (1 not defined because of singularities)
##                                Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                  -0.6287267  0.3082158  -2.040 0.041361 *  
## stage_newcons1                0.3578416  0.5578862   0.641 0.521247    
## stage_newcons2                0.0043790  0.5902286   0.007 0.994080    
## stage_newcons3                0.0925494  0.6807720   0.136 0.891863    
## admission_age                 0.0004099  0.0027326   0.150 0.880768    
## genderM                      -0.0281637  0.0974094  -0.289 0.772484    
## weight_admit                 -0.0276656  0.0036308  -7.620 2.54e-14 ***
## first_stage                   0.2411685  0.1708156   1.412 0.157990    
## stage_newcons1:admission_age  0.0178498  0.0049690   3.592 0.000328 ***
## stage_newcons2:admission_age  0.0219070  0.0053110   4.125 3.71e-05 ***
## stage_newcons3:admission_age  0.0263628  0.0064228   4.105 4.05e-05 ***
## stage_newcons1:genderM        0.4105403  0.1628538   2.521 0.011705 *  
## stage_newcons2:genderM        0.2926366  0.1635858   1.789 0.073633 .  
## stage_newcons3:genderM        0.1316719  0.2070372   0.636 0.524788    
## stage_newcons1:weight_admit  -0.0102460  0.0059355  -1.726 0.084308 .  
## stage_newcons2:weight_admit  -0.0001437  0.0057383  -0.025 0.980016    
## stage_newcons3:weight_admit   0.0100127  0.0068416   1.464 0.143329    
## stage_newcons1:first_stage    0.9043423  0.2256301   4.008 6.12e-05 ***
## stage_newcons2:first_stage    0.7370432  0.2160448   3.412 0.000646 ***
## stage_newcons3:first_stage           NA         NA      NA       NA    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 10571  on 14618  degrees of freedom
## Residual deviance:  8936  on 14600  degrees of freedom
##   (304 observations deleted due to missingness)
## AIC: 8974
## 
## Number of Fisher Scoring iterations: 5

summary(m2_adj1)

## 
## Call:
## glm(formula = mortality_30 ~ stage_newmean * (admission_age + 
##     weight_admit + gender + first_stage), family = binomial, 
##     data = (akis_all_long_complete %>% filter(group == "newmean") %>% 
##         mutate(stage_newmean = as.factor(max_stage))))
## 
## Coefficients: (1 not defined because of singularities)
##                                Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                  -0.1576508  0.3878876  -0.406 0.684424    
## stage_newmean1               -1.0564643  0.6398699  -1.651 0.098726 .  
## stage_newmean2               -0.7006592  0.5468477  -1.281 0.200099    
## stage_newmean3                0.0430433  0.6131998   0.070 0.944039    
## admission_age                -0.0071244  0.0033610  -2.120 0.034030 *  
## weight_admit                 -0.0304899  0.0048119  -6.336 2.35e-10 ***
## genderM                      -0.0625418  0.1255994  -0.498 0.618522    
## first_stage                   0.3992082  0.1300693   3.069 0.002146 ** 
## stage_newmean1:admission_age  0.0205561  0.0056732   3.623 0.000291 ***
## stage_newmean2:admission_age  0.0265655  0.0048815   5.442 5.27e-08 ***
## stage_newmean3:admission_age  0.0262541  0.0055129   4.762 1.91e-06 ***
## stage_newmean1:weight_admit   0.0024194  0.0072056   0.336 0.737044    
## stage_newmean2:weight_admit  -0.0008742  0.0060257  -0.145 0.884645    
## stage_newmean3:weight_admit   0.0076495  0.0065503   1.168 0.242880    
## stage_newmean1:genderM        0.3189261  0.1963144   1.625 0.104255    
## stage_newmean2:genderM        0.2709073  0.1634240   1.658 0.097379 .  
## stage_newmean3:genderM        0.4265736  0.1838586   2.320 0.020335 *  
## stage_newmean1:first_stage    0.9812856  0.1946681   5.041 4.64e-07 ***
## stage_newmean2:first_stage    0.5670661  0.1647149   3.443 0.000576 ***
## stage_newmean3:first_stage           NA         NA      NA       NA    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 10570.6  on 14618  degrees of freedom
## Residual deviance:  8886.5  on 14600  degrees of freedom
##   (304 observations deleted due to missingness)
## AIC: 8924.5
## 
## Number of Fisher Scoring iterations: 6

km_fit_newcons <-
  survfit(Surv(FOLLOWUP_DAYS, DEATH_FLAG) ~ MAX_STAGE_NEW_CONS,
          akis_all_wide)

km_fit_newmean <-
  survfit(Surv(FOLLOWUP_DAYS, DEATH_FLAG) ~ MAX_STAGE_NEW_MEAN,
          akis_all_wide)


aumc_survival_figure_newcons <- km_fit_newcons %>%
ggsurvfit(linewidth = 1) +
  add_confidence_interval() +
  add_quantile() +
  scale_ggsurvfit(x_scales = list(breaks = seq(0, 30, by = 5))) +
  coord_cartesian(xlim = c(0, 30)) +
  theme_classic() +
  scale_fill_discrete(labels=c('0', '1', '2', '3')) +
  scale_color_discrete(labels=c('0', '1', '2', '3')) +
  labs(x="Days", y = "Survival", title = "UOcons",
       color='Maximum KDIGO-UO stage', fill='Maximum KDIGO-UO stage') +
  theme(legend.position = "bottom")

aumc_survival_figure_newmean <- km_fit_newmean %>%
ggsurvfit(linewidth = 1) +
  add_confidence_interval() +
  add_quantile() +
  scale_fill_discrete(labels=c('0', '1', '2', '3')) +
  scale_color_discrete(labels=c('0', '1', '2', '3')) +
  scale_ggsurvfit(x_scales = list(breaks = seq(0, 30, by = 5))) +
  coord_cartesian(xlim = c(0, 30)) +
  theme_classic() +
  labs(x="Days", y = "Survival",title = "UOmean",
       color='Maximum KDIGO-UO stage', fill='Maximum KDIGO-UO stage') +
  theme(legend.position = "bottom")

aumc_kdigo_inter_survival_figure <- ggarrange(
  aumc_survival_figure_newmean,
  aumc_survival_figure_newcons,
  ncol = 1,
  nrow = 2,
  legend = "bottom",
  common.legend = TRUE
) %>% annotate_figure(top = text_grob("AUMCdb", face = "bold", size = 14))

aumc_kdigo_inter_survival_figure

Sensitivity Analysis: Validity Threshold for Durations of Collection

validity_threshold_long <- validity_threshold_wide %>%
  mutate(
         first_stage.cons = FIRST_STAGE_NEW_CONS,
         first_stage.9520 = FIRST_STAGE_NEW_CONS_95_20,
         first_stage.9920 = FIRST_STAGE_NEW_CONS_99_20,
         max_stage.cons = MAX_STAGE_NEW_CONS,
         max_stage.9520 = MAX_STAGE_NEW_CONS_95_20,
         max_stage.9920 = MAX_STAGE_NEW_CONS_99_20,
         .keep = "unused") %>%
  pivot_longer(
    !c(STAY_ID, FOLLOWUP_DAYS, DEATH_FLAG),
    names_sep = "\\.",
    names_to = c(".value", "group")
  ) %>%
  mutate(
    mortality_90 = if_else(FOLLOWUP_DAYS < 91 &
                                 DEATH_FLAG == 1, 1, 0),
    prevalnce_admit = if_else(first_stage > 0, 1, 0),
    Incidence_first_72hr = case_when(first_stage > 0 ~ NA, max_stage == 0 ~ 0, max_stage > 0 ~ 1),
    Incidence_first_72hr_with_stage = ifelse(first_stage == 0 &
                                              max_stage > 0, max_stage, NA),
    .keep = "all"
  ) 

validity_threshold_long %>%
  drop_na(prevalnce_admit) %>%
  transmute(
    group,
    aki_binary = if_else(max_stage > 0, 1, 0),
    max_stage = if_else(max_stage == 0, NA, max_stage),
    prevalnce_admit
  ) %>%
  mutate(group =
           case_when(group == "cons" ~ "No exclusion",
                     group == "9520" ~ "95th precentile for rate bellow 20th precentile",
                     group == "9920" ~ "99th precentile for rate bellow 20th precentile")) %>%
  tbl_summary(
    by = "group",
    missing = "no",
    digits = everything() ~ c(0, 1),
    label = list(
      aki_binary ~ "Oliguric-AKI on the first days",
      prevalnce_admit ~ "Prevalence at admission",
      max_stage ~ "Maximum KDIGO staging"
    )
  )  %>%
  modify_column_indent(columns = label, rows = c(FALSE, TRUE)) %>%
  modify_column_indent(
    columns = label,
    rows = c(FALSE, FALSE, TRUE, TRUE, TRUE),
    double_indent = TRUE
  ) %>%
  add_p()

Characteristic	95th precentile for rate bellow 20th precentile N = 14,8741	99th precentile for rate bellow 20th precentile N = 14,9221	No exclusion N = 14,9321	p-value2
Oliguric-AKI on the first days	4,515 (30.4%)	4,662 (31.2%)	4,690 (31.4%)	0.11
Maximum KDIGO staging				<0.001
1	2,408 (53.3%)	2,458 (52.7%)	2,457 (52.4%)
2	1,658 (36.7%)	1,593 (34.2%)	1,590 (33.9%)
3	449 (9.94%)	611 (13.1%)	643 (13.7%)
Prevalence at admission	969 (6.51%)	1,013 (6.79%)	1,025 (6.86%)	0.4
1 n (%)
2 Pearson’s Chi-squared test

aumc_exclusion_threshold <- validity_threshold_long %>%
  transmute(STAY_ID,
            group,
            first_stage,
            max_stage,
            aki_above_2 = max_stage > 1) %>%
  filter(aki_above_2 == TRUE,
         (group == "cons" | group == "9520" | group == "9920")) %>%
  left_join(table_1, by = "STAY_ID") %>%
      mutate(group =
           case_when(group == "cons" ~ "No exclusion",
                     group == "9520" ~ "95th precentile for rate bellow 20th precentile",
                     group == "9920" ~ "99th precentile for rate bellow 20th precentile")) %>%
  select(
    group,
    admission_age,
    weight_admit,
    gender,
    creat_first,
    creat_peak_72,
    creat_last,
    icu_days,
    rrt_binary,
    hospital_expire_flag
  ) %>%
  tbl_summary(
    by = group,
    type = list(
      c(hospital_expire_flag,
        rrt_binary) ~ "dichotomous",
      c(admission_age,
        weight_admit,
        creat_first,
        creat_peak_72,
        creat_last) ~ "continuous"
    ),
    statistic = c(admission_age,
                  weight_admit,
                  creat_first,
                  creat_peak_72,
                  creat_last) ~ "{mean} ({sd})",
    missing = "no",
    label = list(
      admission_age ~ "Age at Hospital Admission, years",
      gender ~ "Gender",
      weight_admit ~ "Weight at ICU Admission, kg",
      creat_first ~ "First Creatinine in Hospital, mg/dL",
      creat_peak_72 ~ "Peak Creatinine at first days, mg/dL",
      creat_last ~ "Hospital Discharge Creatinine, mg/dL",
      icu_days ~ "Time in ICU, days",
      rrt_binary ~ "Renal replacement therapy",
      hospital_expire_flag ~ "Hospital Mortality"
    )
  ) %>%
  add_p() %>%
  add_stat(
    fns = everything() ~ add_by_n
  ) %>%
  modify_header(starts_with("add_n_stat") ~ "**N**") %>%
  modify_table_body(
    ~ .x %>%
      dplyr::relocate(add_n_stat_1, .before = stat_1) %>%
      dplyr::relocate(add_n_stat_2, .before = stat_2) %>%
      dplyr::relocate(add_n_stat_3, .before = stat_3)
  )

aumc_exclusion_threshold

Characteristic	N	95th precentile for rate bellow 20th precentile N = 2,1071	N	99th precentile for rate bellow 20th precentile N = 2,2041	N	No exclusion N = 2,2331	p-value2
Age at Hospital Admission, years	2,107	67 (14)	2,204	67 (14)	2,233	67 (14)	>0.9
Weight at ICU Admission, kg	2,107	85 (15)	2,204	84 (15)	2,233	84 (15)	>0.9
Gender	2,080		2,175		2,204		>0.9
F		733 (35%)		770 (35%)		776 (35%)
M		1,347 (65%)		1,405 (65%)		1,428 (65%)
First Creatinine in Hospital, mg/dL	2,105	1.67 (1.79)	2,202	1.71 (1.85)	2,231	1.72 (1.87)	0.8
Peak Creatinine at first days, mg/dL	2,095	2.43 (1.86)	2,191	2.45 (1.90)	2,219	2.45 (1.91)	>0.9
Hospital Discharge Creatinine, mg/dL	2,105	1.66 (1.34)	2,202	1.67 (1.36)	2,231	1.67 (1.35)	>0.9
Time in ICU, days	2,107	5 (3, 13)	2,204	5 (2, 13)	2,233	5 (2, 13)	>0.9
Renal replacement therapy	2,107	646 (31%)	2,204	687 (31%)	2,233	696 (31%)	>0.9
Hospital Mortality	2,107	654 (31%)	2,204	689 (31%)	2,233	698 (31%)	>0.9
1 Mean (SD); n (%); Median (Q1, Q3)
2 Kruskal-Wallis rank sum test; Pearson’s Chi-squared test

---

Clinical Outcomes

Describing the prevalence of oliguric-AKI upon admission and incidence at the first ICU day

akis_all_long %>%
  filter(group == "newcons") %>%
  select(prevalnce_admit,
         Incidence_first_72hr,
         max_stage) %>%
  drop_na(prevalnce_admit) %>%
  transmute(
    aki_binary = if_else(max_stage > 0, 1, 0),
    max_stage = if_else(max_stage == 0, NA, max_stage),
    prevalnce_admit
  ) %>%
  tbl_summary(
    missing = "no",
    digits = everything() ~ c(0, 1),
    label = list(
      aki_binary ~ "Oliguric-AKI on the first days",
      prevalnce_admit ~ "Prevalence at admission",
      max_stage ~ "Maximum KDIGO staging"
    )
  ) %>%
  modify_column_indent(columns = label, 
                       rows = c(FALSE, TRUE)) %>%
  modify_column_indent(columns = label, 
                       rows = c(FALSE, FALSE, TRUE, TRUE, TRUE),
                       double_indent = TRUE)

Characteristic	N = 14,9231
Oliguric-AKI on the first days	4,688 (31.4%)
Maximum KDIGO staging
1	2,456 (52.4%)
2	1,589 (33.9%)
3	643 (13.7%)
Prevalence at admission	1,024 (6.86%)
1 n (%)

aki_uo_analysis <- left_join(akis_all_wide, uo_ml_kg_hr, by = "STAY_ID") %>%
  drop_na(FIRST_POSITIVE_STAGE_UO_CONS_TIME, 
          TIME_INTERVAL_FINISH, 
          MAX_STAGE_NEW_CONS) %>%
  transmute(STAY_ID,
         MAX_STAGE_NEW_CONS = as.character(MAX_STAGE_NEW_CONS),
         FIRST_POSITIVE_STAGE_UO_CONS_TIME,
         TIME_INTERVAL_FINISH,
         UO_KG = ML_KG_HR
         ) %>%
  mutate(TIME = as.double(difftime(TIME_INTERVAL_FINISH, 
                                   FIRST_POSITIVE_STAGE_UO_CONS_TIME, 
                                   units = c("hour")))) %>%
  filter(TIME >= -48 & TIME <= 48)

First Oliguric-AKI Events

table 1 for ICU stays with identified oliguric AKI in the first 72 hours of admission, stratified by max kdigo-uo stage (ICU stays with AKI at admission were excluded):

table1_akis <- akis_all_wide  %>%
  filter(MAX_STAGE_NEW_CONS >= 0) %>%
  select(STAY_ID, MAX_STAGE_NEW_CONS) %>%
  left_join(table_1, by = "STAY_ID")

table_1_staging <- table1_akis %>%
  select(
    MAX_STAGE_NEW_CONS,
    admission_age,
    weight_admit,
    gender,
    creat_first,
    creat_peak_72,
    creat_last,
    icu_days,
    rrt_binary,
    hospital_expire_flag
  ) %>%
  mutate(
    staging = case_when(
      MAX_STAGE_NEW_CONS == 0 ~ "No AKI",
      MAX_STAGE_NEW_CONS == 1 ~ "Stage 1",
      MAX_STAGE_NEW_CONS == 2 ~ "Stage 2",
      MAX_STAGE_NEW_CONS == 3 ~ "Stage 3",
      .default = NA
    ),
    .keep = "unused"
  ) %>%
  tbl_summary(
    by = staging,
    type = list(
      c(hospital_expire_flag, rrt_binary) ~ "dichotomous",
      c(
        admission_age,
        weight_admit,
        creat_peak_72,
        creat_first,
        creat_last
      ) ~ "continuous"
    ),
    statistic = c(
      admission_age,
      weight_admit,
      creat_peak_72,
      creat_first,
      creat_last
    ) ~ "{mean} ({sd})",
    digits = icu_days ~ c(1,1),
    missing = "no",
    missing_text = "-",
    label = list(
      admission_age ~ "Age at Hospital Admission, years",
      gender ~ "Gender",
      weight_admit ~ "Weight at ICU Admission, kg",
      creat_first ~ "First Creatinine in Hospital, mg/dL",
      creat_peak_72 ~ "Peak Creatinine at first days, mg/dL",
      creat_last ~ "Hospital Discharge Creatinine, mg/dL",
      icu_days ~ "Time in ICU, days",
      rrt_binary ~ "Renal replacement therapy",
      hospital_expire_flag ~ "Hospital Mortality"
    )
  ) %>%
  add_p(test = c(
    admission_age,
    weight_admit,
    creat_peak_72,
    creat_first,
    creat_last
  ) ~ "aov")

table_1_staging

Characteristic	No AKI N = 10,2351	Stage 1 N = 2,4561	Stage 2 N = 1,5891	Stage 3 N = 6431	p-value2
Age at Hospital Admission, years	62 (15)	65 (15)	68 (14)	66 (15)	<0.001
Weight at ICU Admission, kg	79 (14)	83 (15)	85 (15)	83 (16)	<0.001
Gender					<0.001
F	3,064 (31%)	879 (36%)	544 (35%)	231 (37%)
M	6,933 (69%)	1,540 (64%)	1,028 (65%)	400 (63%)
First Creatinine in Hospital, mg/dL	1.07 (0.65)	1.19 (0.93)	1.44 (1.50)	2.40 (2.42)	<0.001
Peak Creatinine at first days, mg/dL	1.06 (0.64)	1.29 (0.94)	1.92 (1.59)	3.74 (1.98)	<0.001
Hospital Discharge Creatinine, mg/dL	0.95 (0.50)	1.07 (0.81)	1.41 (1.20)	2.30 (1.50)	<0.001
Time in ICU, days	1.0 (0.9, 2.0)	2.5 (1.0, 6.9)	4.7 (2.1, 11.3)	7.2 (3.0, 18.2)	<0.001
Renal replacement therapy	44 (0.4%)	61 (2.5%)	213 (13%)	483 (75%)	<0.001
Hospital Mortality	424 (4.1%)	310 (13%)	389 (24%)	308 (48%)	<0.001
1 Mean (SD); n (%); Median (Q1, Q3)
2 One-way analysis of means; Pearson’s Chi-squared test; Kruskal-Wallis rank sum test

UO onset at UOcons event:

aumc_uo_cons_figure <- aki_uo_analysis %>% 
ggplot(aes(TIME, UO_KG, color=MAX_STAGE_NEW_CONS, fill=MAX_STAGE_NEW_CONS))  + 
           # linetype=MAX_STAGE_NEW_CONS))  + 
  geom_hline(yintercept=0.3, size = 0.3, color = "#cccccc") +
  geom_hline(yintercept=0.5, size = 0.3, color = "#cccccc") +
  geom_vline(xintercept=0, size = 0.3, color = "black", linetype = "dashed") +
  stat_summary(fun = median, geom="line") +
  scale_x_continuous(breaks = seq(-24, 48, by=6)) +
  scale_y_continuous(breaks = c(0, 0.3, 0.5)) +
  coord_cartesian(xlim = c(-12, 24), ylim = c(0, 1.7)) +
  # xlim(-24, 48) +
  stat_summary(fun.min = function(z) { quantile(z,0.25) },
               fun.max = function(z) { quantile(z,0.75) },
               geom="ribbon", colour = NA, alpha=0.2) +
  labs(x="Time around AKI onset (hour)", y = "Urine output (ml/kg/hr)", 
       color="Maximum KDIGO-UO stage", fill="Maximum KDIGO-UO stage") + 
  theme_classic() + # remove panel background and gridlines
  scale_color_manual(values = pal_jama("default")(4)[2:4]) +
  scale_fill_manual(values = pal_jama("default")(4)[2:4]) +
  theme(
    legend.position = "none"
  )

aumc_uo_cons_figure

 

Survival Analysis

km_fit <- survfit2(Surv(FOLLOWUP_DAYS, mortality_30) ~ MAX_STAGE_NEW_CONS, data = akis_all_wide)

aumc_survival_cons_figure <- km_fit %>%
ggsurvfit(linewidth = 1) +
  scale_ggsurvfit(x_scales = list(breaks = seq(0, 30, by = 5))) +
  coord_cartesian(xlim = c(-1, 31)) +
  theme_classic() +
  labs(x="Days", y = "Survival", 
       color='Maximum KDIGO-UO stage', fill='Maximum KDIGO-UO stage') +
  scale_color_jama() +
  scale_fill_jama() +
  theme(legend.position = "bottom",
        plot.title = element_text(color = "#0099F8", size = 16, face = "bold"),
        plot.subtitle = element_text(size = 10, face = "bold")) +
  add_risktable() +
  add_pvalue(caption = "Log-rank {p.value}")

aumc_survival_cons_figure

table of survival probabilities:

aumc_survival_table <-
  km_fit %>% tbl_survfit(times = c(7, 30, 90, 365),
                                                 label = "Maximum KDIGO-UO stage",
                                                 label_header = "**Day {time}**")
aumc_survival_table

Characteristic	Day 7	Day 30	Day 90	Day 365
Maximum KDIGO-UO stage
0	97% (96%, 97%)	94% (94%, 95%)	94% (94%, 95%)	94% (94%, 95%)
1	90% (89%, 91%)	84% (83%, 86%)	84% (83%, 86%)	84% (83%, 86%)
2	82% (80%, 84%)	73% (71%, 75%)	73% (71%, 75%)	73% (71%, 75%)
3	66% (62%, 70%)	51% (48%, 55%)	51% (48%, 55%)	51% (48%, 55%)

Log rank for each pair:

survdiff(Surv(FOLLOWUP_DAYS, mortality_30) ~ MAX_STAGE_NEW_CONS, data = akis_all_wide)

## Call:
## survdiff(formula = Surv(FOLLOWUP_DAYS, mortality_30) ~ MAX_STAGE_NEW_CONS, 
##     data = akis_all_wide)
## 
## n=14923, 381 observations deleted due to missingness.
## 
##                          N Observed Expected (O-E)^2/E (O-E)^2/V
## MAX_STAGE_NEW_CONS=0 10235      605   1221.6     311.2    1070.7
## MAX_STAGE_NEW_CONS=1  2456      384    277.9      40.5      48.7
## MAX_STAGE_NEW_CONS=2  1589      428    169.3     395.1     442.0
## MAX_STAGE_NEW_CONS=3   643      312     60.2    1053.9    1103.5
## 
##  Chisq= 1820  on 3 degrees of freedom, p= <2e-16

akis_all_wide_non_01 <- akis_all_wide %>%
  filter(MAX_STAGE_NEW_CONS == 0 | MAX_STAGE_NEW_CONS == 1)
survdiff(Surv(FOLLOWUP_DAYS, mortality_30) ~ MAX_STAGE_NEW_CONS, data = akis_all_wide_non_01)

## Call:
## survdiff(formula = Surv(FOLLOWUP_DAYS, mortality_30) ~ MAX_STAGE_NEW_CONS, 
##     data = akis_all_wide_non_01)
## 
##                          N Observed Expected (O-E)^2/E (O-E)^2/V
## MAX_STAGE_NEW_CONS=0 10235      605      806        50       271
## MAX_STAGE_NEW_CONS=1  2456      384      183       220       271
## 
##  Chisq= 272  on 1 degrees of freedom, p= <2e-16

akis_all_wide_non_12 <- akis_all_wide %>%
  filter(MAX_STAGE_NEW_CONS == 1 | MAX_STAGE_NEW_CONS == 2)
survdiff(Surv(FOLLOWUP_DAYS, mortality_30) ~ MAX_STAGE_NEW_CONS, data = akis_all_wide_non_12)

## Call:
## survdiff(formula = Surv(FOLLOWUP_DAYS, mortality_30) ~ MAX_STAGE_NEW_CONS, 
##     data = akis_all_wide_non_12)
## 
##                         N Observed Expected (O-E)^2/E (O-E)^2/V
## MAX_STAGE_NEW_CONS=1 2456      384      504      28.7        77
## MAX_STAGE_NEW_CONS=2 1589      428      308      47.1        77
## 
##  Chisq= 77  on 1 degrees of freedom, p= <2e-16

akis_all_wide_non_23 <- akis_all_wide %>%
  filter(MAX_STAGE_NEW_CONS == 2 | MAX_STAGE_NEW_CONS == 3)
survdiff(Surv(FOLLOWUP_DAYS, mortality_30) ~ MAX_STAGE_NEW_CONS, data = akis_all_wide_non_23)

## Call:
## survdiff(formula = Surv(FOLLOWUP_DAYS, mortality_30) ~ MAX_STAGE_NEW_CONS, 
##     data = akis_all_wide_non_23)
## 
##                         N Observed Expected (O-E)^2/E (O-E)^2/V
## MAX_STAGE_NEW_CONS=2 1589      428      546      25.6       101
## MAX_STAGE_NEW_CONS=3  643      312      194      72.2       101
## 
##  Chisq= 101  on 1 degrees of freedom, p= <2e-16

# aumc_aki_epi <- aki_epi
# aumc_akis <- akis
aumc_t1a <- t1a
aumc_t1b <- t1b
aumc_table_1 <- table_1
aumc_table_1_akis <- table1_akis
aumc_uo_rate <- uo_rate
aumc_aki_epi <- akis_all_long %>%
  filter(group == "newcons")
aumc_akis_all_wide <- akis_all_wide
aumc_table_1_staging <- table_1_staging
aumc_rr_table_data <- rr_table_data
aumc_rr_table_data_adj <- rr_table_data %>% 
  left_join(bind_rows(rr_table_data_m1_adj1, rr_table_data_m2_adj1))

save(aumc_t1a,
     aumc_t1b,
     aumc_table_1,
     aumc_table_1_akis,
     aumc_uo_rate,
     aumc_uo_cons_figure,
     aumc_survival_cons_figure,
     aumc_survival_table,
     aumc_aki_epi,
     aumc_akis_all_wide,
     aumc_table_1_staging,
     aumc_rr_table_data,
     aumc_rr_table_data_adj,
     file = "paper_aumc.Rda")

aumc_uo_rate <- uo_rate
aumc_hourly_uo <- hourly_uo
aumc_raw_uo_eligible <- raw_uo_eligible

save(
  bS3a,
  bS6_a,
  bS6_b,
  # bS7_a,
  # bS7_b,
  bS7_c,
  bS7_d,
  # bS8_a,
  # bS8_b,
  bS8_c,
  bS8_d,
  bS8_e,
  # bS8_f,
  bS9_a,
  bS9_b,
  bS9_c,
  bS9_d,
  bS11,
  # aumc_uo_rate,
  # aumc_hourly_uo,
  # aumc_raw_uo_eligible,
  aumc_exclusion_threshold,
  aumc_kdigo_inter_aki_table,
  aumc_kdigo_inter_cons_mean,
  aumc_kdigo_inter_bic,
  aumc_kdigo_inter_survival_table,
  aumc_kdigo_inter_survival_table_adj,
  aumc_kdigo_inter_survival_figure,
  file = "bs_data.Rda"
)

---

Technical Details

R Session Info:

sessionInfo()

## R version 4.4.1 (2024-06-14)
## Platform: aarch64-apple-darwin20
## Running under: macOS Sonoma 14.5
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRblas.0.dylib 
## LAPACK: /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRlapack.dylib;  LAPACK version 3.12.0
## 
## 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: Asia/Jerusalem
## tzcode source: internal
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] marginaleffects_0.21.0 rms_6.8-1              Hmisc_5.1-3           
##  [4] gt_0.11.0              ggsurvfit_1.1.0        ggsci_3.2.0           
##  [7] gtsummary_2.0.0        nortest_1.0-4          survminer_0.4.9       
## [10] ggpubr_0.6.0           survival_3.7-0         rmdformats_1.0.4      
## [13] kableExtra_1.4.0       broom_1.0.6            quantreg_5.98         
## [16] SparseM_1.84-2         rlang_1.1.4            ggforce_0.4.2         
## [19] ggpmisc_0.6.0          ggpp_0.5.8-1           scales_1.3.0          
## [22] ggbreak_0.1.2          psych_2.4.6.26         finalfit_1.0.8        
## [25] reshape2_1.4.4         lubridate_1.9.3        forcats_1.0.0         
## [28] stringr_1.5.1          dplyr_1.1.4            purrr_1.0.2           
## [31] readr_2.1.5            tidyr_1.3.1            tibble_3.2.1          
## [34] ggplot2_3.5.1          tidyverse_2.0.0        bigrquery_1.5.1       
## [37] DBI_1.2.3              pacman_0.5.1          
## 
## loaded via a namespace (and not attached):
##   [1] splines_4.4.1        polspline_1.1.25     ggplotify_0.1.2     
##   [4] polyclip_1.10-7      rpart_4.1.23         lifecycle_1.0.4     
##   [7] rstatix_0.7.2        lattice_0.22-6       MASS_7.3-61         
##  [10] insight_0.20.2       backports_1.5.0      magrittr_2.0.3      
##  [13] sass_0.4.9           rmarkdown_2.27       jquerylib_0.1.4     
##  [16] yaml_2.3.10          askpass_1.2.0        cowplot_1.1.3       
##  [19] RColorBrewer_1.1-3   minqa_1.2.7          multcomp_1.4-26     
##  [22] abind_1.4-5          clock_0.7.1          yulab.utils_0.1.5   
##  [25] nnet_7.3-19          TH.data_1.1-2        tweenr_2.0.3        
##  [28] rappdirs_0.3.3       sandwich_3.1-0       labelled_2.13.0     
##  [31] KMsurv_0.1-5         cards_0.2.0          MatrixModels_0.5-3  
##  [34] cardx_0.2.0          svglite_2.1.3        commonmark_1.9.1    
##  [37] codetools_0.2-20     xml2_1.3.6           tidyselect_1.2.1    
##  [40] shape_1.4.6.1        aplot_0.2.3          farver_2.1.2        
##  [43] lme4_1.1-35.5        base64enc_0.1-3      broom.helpers_1.15.0
##  [46] jsonlite_1.8.8       mitml_0.4-5          Formula_1.2-5       
##  [49] iterators_1.0.14     systemfonts_1.1.0    foreach_1.5.2       
##  [52] tools_4.4.1          Rcpp_1.0.13          glue_1.7.0          
##  [55] mnormt_2.1.1         gridExtra_2.3        pan_1.9             
##  [58] mgcv_1.9-1           xfun_0.46            withr_3.0.0         
##  [61] fastmap_1.2.0        boot_1.3-30          fansi_1.0.6         
##  [64] openssl_2.2.0        digest_0.6.36        timechange_0.3.0    
##  [67] R6_2.5.1             gridGraphics_0.5-1   mice_3.16.0         
##  [70] colorspace_2.1-1     markdown_1.13        utf8_1.2.4          
##  [73] generics_0.1.3       data.table_1.15.4    httr_1.4.7          
##  [76] htmlwidgets_1.6.4    pkgconfig_2.0.3      gtable_0.3.5        
##  [79] survMisc_0.5.6       brio_1.1.5           htmltools_0.5.8.1   
##  [82] carData_3.0-5        bookdown_0.40        png_0.1-8           
##  [85] ggfun_0.1.5          knitr_1.48           km.ci_0.5-6         
##  [88] rstudioapi_0.16.0    tzdb_0.4.0           checkmate_2.3.1     
##  [91] nlme_3.1-165         curl_5.2.1           nloptr_2.1.1        
##  [94] cachem_1.1.0         zoo_1.8-12           parallel_4.4.1      
##  [97] foreign_0.8-87       pillar_1.9.0         grid_4.4.1          
## [100] vctrs_0.6.5          car_3.1-2            jomo_2.7-6          
## [103] xtable_1.8-4         cluster_2.1.6        htmlTable_2.4.3     
## [106] evaluate_0.24.0      mvtnorm_1.2-5        cli_3.6.3           
## [109] compiler_4.4.1       ggsignif_0.6.4       labeling_0.4.3      
## [112] plyr_1.8.9           fs_1.6.4             stringi_1.8.4       
## [115] viridisLite_0.4.2    munsell_0.5.1        glmnet_4.1-8        
## [118] Matrix_1.7-0         hms_1.1.3            patchwork_1.2.0     
## [121] bit64_4.0.5          haven_2.5.4          highr_0.11          
## [124] gargle_1.5.2         memoise_2.0.1        bslib_0.7.0         
## [127] bit_4.0.5            polynom_1.4-1