In windyzn/urinaryDBP:

# Only set if the Rmd file is not in the parent directory (ie. 'projectname/')
knitr::opts_knit$set(root.dir = '../')

knitr::opts_chunk$set(collapse = TRUE, echo = FALSE, message = FALSE, warning = FALSE)
library(tidyverse)
library(PROMISE.data)
devtools::load_all()
# load_data(update = TRUE)
set_options()
# extrafont::loadfonts(device="win")

# source('.Rprofile')
# run_setup()
# load_data('data/project_data.rda')
# dim(ds)
ds <- PROMISE.data::PROMISE%>%

  PROMISE.data::combine_datasets(ogtt)%>%

  # Data wrangling
    dplyr::mutate(
      UDBP = ifelse(UDBP == 0.01, NA, UDBP),
      UrineCreatinine = ifelse(SID == 2028, 9, UrineCreatinine),
      ACR = round(UrineMicroalbumin / UrineCreatinine, digits = 2),
      Ethnicity = as.character(Ethnicity),
      isAfrican = ifelse(Ethnicity == 'African', 1, 0),
      Ethnicity = ifelse(
        Ethnicity %in% c('African', 'First Nations', 'Other'),
        'Other',
        Ethnicity
      ),
      fVN = factor(
        VN,
        levels = c(0, 3, 6),
        labels = c("Baseline", "3Year", "6Year"),
        ordered = TRUE
      ),
      YearsFromBaseline = MonthsFromBaseline/12,
      fMedsBP = factor(
        MedsBloodPressure,
        levels = c(0, 1),
        labels = c("No", "Yes")
      ),
      Hypertension = ifelse(Systolic > 140 | Diastolic > 90, 1,
                                   0),
      Hypertension = factor(Hypertension,
                            levels = c(0, 1),
                            labels = c("No", "Yes")),
      dmStatus = ifelse(DM == 1, "DM",
                         ifelse(IFG == 1 |
                                  IGT == 1, "PreDM",
                                "NGT")),
      dmStatus = ifelse(is.na(dmStatus), "NGT", dmStatus),
      dmStatus = factor(dmStatus,
                        levels = c("NGT", "PreDM", "DM"),
                        ordered = TRUE),
      acrStatus = ifelse(
        ACR < 2,
        'Normoalbuminuria',
        ifelse(ACR > 20, 'Macroalbuminuria',
               "Microalbuminuria")
      ),
      acrStatus = factor(
        acrStatus,
        levels = c("Normoalbuminuria", "Microalbuminuria", "Macroalbuminuria"),
        ordered = TRUE
      ),
      creat.mgdl = Creatinine * 0.011312,
      nSex = ifelse(Sex == "Female", 0, 1),
      eGFR = nephro::CKDEpi.creat(creat.mgdl, nSex, Age, isAfrican),
      eGFRStatus = cut(eGFR,
                       breaks = c(-Inf, 60, 90, Inf),
                       labels = c("Moderate", "Mild", "Normal"),
                       ordered_result = TRUE),
      # eGFRStatus = factor(eGFRStatus,
      #                     levels = c("Normal", "Mild", "Moderate"),
      #                     ordered = TRUE),
      udbpStatus = cut(UDBP,
                       breaks = c(0, 1.23, 60, Inf),
                       labels = c("Trace", "Normal", "High"),
                       ordered_result = TRUE),
      udbpTertile = ntile(UDBP, 3),
      udbpCr = UDBP / UrineCreatinine,
      logudbpCr = log(udbpCr),
      udbpCrTertile = ntile(udbpCr, 3),
      vitdStatus = cut(VitaminD,
                       breaks = c(-Inf, 50, 75, Inf),
                       labels = c("Deficient", "Insufficient", "Sufficient"),
                       ordered = TRUE),
      Season = ifelse(lubridate::month(VisitDate) %in% c("5", "6", "7", "8", "9", "10"),
                      "Summer", "Winter")
    ) %>%
    dplyr::filter(UDBP < 10000) %>%
    dplyr::filter(eGFR < 200) %>%
    dplyr::select(
      SID,
      VN,
      fVN,
      MonthsFromBaseline,
      YearsFromBaseline,
      VisitDate,
      BMI,
      Waist,
      Age,
      Sex,
      Ethnicity,
      Glucose0_OGTT,
      Glucose120_OGTT,
      Dysgly,
      DM,
      IFG,
      IGT,
      dmStatus,
      acrStatus,
      eGFRStatus,
      udbpStatus,
      udbpTertile,
      udbpCrTertile,
      eGFR,
      ACR,
      UrineMicroalbumin,
      UrineCreatinine,
      Creatinine,
      UDBP,
      udbpCr,
      logudbpCr,
      VitaminD,
      vitdStatus,
      Season,
      MeanArtPressure,
      Systolic,
      Diastolic,
      Hypertension,
      PTH,
      ALT,
      fMedsBP,
      MET,
      CRP,
      Canoe,
      MedsBloodPressure
      # dplyr::matches("meds")
    )

dsBase <- ds %>% 
  dplyr::filter(fVN == "Baseline")

# Subjects with measurements at all visit numbers

dsComplete <- ds %>%
  dplyr::group_by(SID) %>%
  dplyr::filter(n() == 3) %>%
  dplyr::ungroup()

## Include captions below using `captioner` package

fig_nums <- captioner::captioner(prefix = 'FIGURE')
supfig_nums <- captioner::captioner(prefix = 'SUPPLEMENTARY FIGURE')
tab_nums <- captioner::captioner(prefix = 'TABLE')
suptab_nums <- captioner::captioner(prefix = 'SUPPLEMENTARY TABLE')

PAPER 1: KIDNEY DYSFUNCTION RESULTS

Subject Characteristics

# tab_nums(kidney_subchar_baseline, "Subject characteristics across uVDBP tertiles at baseline")

dsBaseNoBPMeds <- dsBase %>% 
  dplyr::filter(fMedsBP == "No")

subchar_table_tert <- tableone::CreateTableOne(
  vars = c("Age",
           "Sex",
           "Ethnicity",
           "BMI",
           "Waist",
           "eGFR",
           "ACR",
           "UrineCreatinine",
           "UrineMicroalbumin",
           "UrinaryCalcium",
           "UDBP",
           "udbpCr",
           "Systolic",
           "Diastolic",
           "MeanArtPressure",
           "dmStatus"),
  strata = c("udbpTertile"),
  data = dsBaseNoBPMeds,
  factorVars = c("Sex", "Ethnicity", "dmStatus")
) %>% 
  print(nonnormal = c("UDBP",
                      "ACR",
                      "UrineMicroalbumin"),
        quote = FALSE,
        noSpaces = TRUE)

# explore non-linear variables (delete later)

dsBase %>% 
  scatter_plot("udbpCr", "BMI",
               "uVDBP:cr", "BMI")

# explore (delete later)

ds %>% 
  dplyr::select(VN, SID, UDBP) %>% 
  dplyr::filter(UDBP < 1.23) %>% 
  summarise(n = n())

TABLE 2: Subject characteristics across UDBP:cR tertiles

subchar_table_tert <- tableone::CreateTableOne(
  vars = c("Age",
           "Sex",
           "Ethnicity",
           "BMI",
           "Waist",
           "eGFR",
           "ACR",
           "UrineCreatinine",
           "UrineMicroalbumin",
           "UrinaryCalcium",
           "UDBP",
           "udbpCr",
           "Systolic",
           "Diastolic",
           "MeanArtPressure",
           "dmStatus"),
  # strata = c("udbpCrTertile"),
  data = dsBase,
  factorVars = c("Sex", "Ethnicity", "dmStatus")
) %>% 
  print(nonnormal = c("UDBP",
                      "ACR",
                      "UrineMicroalbumin"),
        quote = FALSE,
        noSpaces = TRUE)

TABLE 3: Subject characteristics across time

subchar_table_time <- tableone::CreateTableOne(
  vars = c("Age",
           "Sex",
           "Ethnicity",
           "BMI",
           "Waist",
           "eGFR",
           "ACR",
           "UrineCreatinine",
           "UrineMicroalbumin",
           "UDBP",
           "udbpCr",
           "Systolic",
           "Diastolic",
           "MeanArtPressure",
           "dmStatus"),
  strata = c("fVN"),
  data = ds,
  factorVars = c("Sex", "Ethnicity", "dmStatus")
) %>% 
  print(nonnormal = c("UDBP",
                      "udbpCr",
                      "ACR",
                      "UrineMicroalbumin"),
        quote = FALSE,
        noSpaces = TRUE)

TABLE 4: Subject characteristics across time (complete data)

subchar_table_time <- tableone::CreateTableOne(
  vars = c("Age",
           "Sex",
           "Ethnicity",
           "BMI",
           "Waist",
           "eGFR",
           "ACR",
           "UrineCreatinine",
           "UrineMicroalbumin",
           "UDBP",
           "udbpCr",
           "Systolic",
           "Diastolic",
           "MeanArtPressure",
           "dmStatus"),
  strata = c("fVN"),
  data = dsComplete,
  factorVars = c("Sex", "Ethnicity", "dmStatus")
) %>% 
  print(nonnormal = c("UDBP",
                      "udbpCr",
                      "ACR",
                      "UrineMicroalbumin"),
        quote = FALSE,
        noSpaces = TRUE)

Incidence

# desc_incidence() works for 0, 1 variables, not factors

ds %>% 
  dplyr::mutate(microalbumin = ifelse(acrStatus == "Microalbuminuria", 1, 0),
                macroalbumin = ifelse(acrStatus == "Macroalbuminuria", 1, 0)) %>% 
  dplyr::select(SID, VN, DM, microalbumin, macroalbumin) %>% 
  PROMISE.desc::desc_incidence("microalbumin")

Cross-sectional

ACR

# Clean data
acr <- ds %>% 
  dplyr::filter(fVN == "Baseline") %>% 
  dplyr::mutate(acrStatus2 = ifelse(acrStatus == "Normoalbuminuria", "Normoalbuminuria",
                                    "Albuminuria"),
                acrStatus2 = factor(acrStatus2,
                                    levels = c("Normoalbuminuria", "Albuminuria"),
                                    ordered = TRUE),
                acrLabels = ifelse(acrStatus == "Normoalbuminuria", "Normoalbuminuria \n(<2 mg/mmol)",
                                   ifelse(acrStatus == "Microalbuminuria", "Microalbuminuria \n(2-20 mg/mmol)",
                                          "Macroalbuminuria \n(>20 mg/mmol)")),
                acrLabels = factor(acrLabels,
                                   levels = c("Normoalbuminuria \n(<2 mg/mmol)",
                                              "Microalbuminuria \n(2-20 mg/mmol)",
                                              "Macroalbuminuria \n(>20 mg/mmol)"),
                                   ordered = TRUE)) %>% 
  dplyr::select(acrStatus, acrStatus2, acrLabels, UDBP, udbpCr) %>%
  na.omit()

# Box plot of uVDBP in different albuminuria categories
acr %>% 
  box_plot("acrStatus2", "log(udbpCr)", 
            "Albuminuria",
            "log uVDBP:cr (μg/mmol)")

# n and values
acr %>% 
  group_by(acrStatus2) %>% 
  summarise(n = n(),
            median = median(log(udbpCr)),
            lower = quantile(log(udbpCr), probs=0.25),
            higher = quantile(log(udbpCr), probs=0.75))


give.n <- function(x){ return(c(ymin = median, label = length(x))) }
boxplot = boxplot + stat_summary(fun.data = give.n, geom = "text", fun.y = median)

# ANOVA
anova <- aov(formula = log(udbpCr)~acrStatus2, data = acr)
summary(anova)
TukeyHSD(anova)
rm(anova)

# Scatterplot of ACR and uVDBP ----------------------------------

dsBase %>% 
  scatter_plot_transparent("log(ACR)", "log(udbpCr)",
               "log albumin:cr ratio (mg/mmol)",
               "log uVDBP:cr (μg/mmol)")

# Spearman Correlation ------------------------------------------

ds %>% 
  filter(fVN == "Baseline") %>% 
  cor.test(formula = ~ ACR + udbpCr, data = ., method = "spearman")

# Linear Regression ---------------------------------------------

dsBase %>% 
  prep_mason_data() %>% 
  mason_glm(y = "ACR",
            x = "udbpCr",
            covars = c("Age", "Sex", "Ethnicity", "dmStatus")
            ) %>% 
  gee_results_table() %>% 
  pander::pander()

eGFR

# Hyperfiltration (eGFR >124.46 at baseline, >121.80 at 3-year, >120.14 at 6-year)
ds %>% 
  dplyr::filter(fVN == "6Year") %>% 
  # dplyr::summarise(mean = mean(eGFR),
  #                  sd = sd(eGFR))
  dplyr::filter(eGFR > 120.14) %>% 
  summarise(n = n())

# Clean data
eGFR <- dsBase %>% 
  dplyr::mutate(eGFRLabels = ifelse(eGFRStatus == "Normal", "Normal \n(>90 mL/min/1.73m²)",
                                   ifelse(eGFRStatus == "Mild", "Mild \n(60-89 mL/min/1.73m²)",
                                          "Moderate \n(<60 mL/min/1.73m²)")),
                eGFRLabels = factor(eGFRLabels,
                                   levels = c("Normal \n(>90 mL/min/1.73m²)",
                                              "Mild \n(60-89 mL/min/1.73m²)",
                                              "Moderate \n(<60 mL/min/1.73m²)"),
                                   ordered = TRUE)) %>%
  select(eGFRStatus, eGFRLabels, udbpCr, UDBP) %>%
  na.omit()

# Box plot of uVDBP in different eGFR categories
eGFR %>% 
  box_plot_slides("eGFRLabels", "log(udbpCr)", 
            "Kidney Status",
            "log uVDBP:cr (μg/mmol)")

# n and values
eGFR %>% 
  group_by(eGFRStatus) %>% 
  summarise(n = n(),
            median = median(log(udbpCr)),
            lower = quantile(log(udbpCr), probs=0.25),
            higher = quantile(log(udbpCr), probs=0.75))

# ANOVA
anova <- aov(formula = log(UDBP)~eGFRStatus, data = eGFR)
summary(anova)
TukeyHSD(anova)
rm(anova)

# Scatterplot of eGFR and uVDBP ----------------------------------

dsBase %>% 
  scatter_plot_transparent("log(eGFR)", "log(udbpCr)",
               "log eGFR (ml/min/1.73m^2)",
               "log uVDBP:cr (μg/mmol)")

# Spearman Correlation ------------------------------------------

ds %>% 
  filter(fVN == "Baseline") %>% 
  cor.test(formula = ~ eGFR + udbpCr, data = ., method = "spearman")

# Linear Regression ---------------------------------------------

dsBase %>% 
  prep_mason_data() %>% 
  mason_glm(y = "eGFR",
            x = "UDBP",
            covars = c("Age", "Sex", "Ethnicity", "dmStatus")
            ) %>% 
  gee_results_table(table = FALSE)

ds %>% 
  prep_mason_data() %>% 
  mason_glm(y = "eGFR",
            x = "udbpCr"
            # covars = "ageBase"
            )

uVDBP vs uVDBP:cr

# Scatterplot of uVDBP and uVDBP:cr -----------------------------

ds %>% 
  scatter_plot("log(UDBP)", "log(udbpCr)",
               "log uVDBP",
               "log uVDBP:cr", facet = TRUE)

ds %>% 
  scatter_plot("log(UDBP)", "log(UrineCreatinine)",
               "log uVDBP",
               "log uCreatinine", facet = TRUE)

# Spearman Correlation ------------------------------------------

ds %>% 
  filter(fVN == "Baseline") %>% 
  cor.test(formula = ~ udbpCr + UDBP, data = ., method = "spearman")

# Scatterplot of eGFR and uVDBP:cr -----------------------------

ds %>% 
  dplyr::filter(fVN == "Baseline") %>% 
  scatter_plot("log(eGFR)", "log(UDBP)",
               "log eGFR",
               "log uVDBP")

# Spearman Correlation ------------------------------------------

ds %>% 
  filter(fVN == "Baseline") %>% 
  cor.test(formula = ~ UDBP + eGFR, data = ., method = "spearman")

# Scatterplot of eGFR and uCr -----------------------------

ds %>% 
  dplyr::filter(fVN == "6Year") %>% 
  scatter_plot("log(eGFR)", "log(UrineCreatinine)",
               "log eGFR",
               "log creatinine")

# Spearman Correlation ------------------------------------------

ds %>% 
  filter(fVN == "6Year") %>% 
  cor.test(formula = ~ UrineCreatinine + eGFR, data = ., method = "spearman")

Medication

There are 705 values across all time points with blood pressure medication data
This is approximately half of all observations (1852 vs 705)
CHECK IF MISSINGNESS == NOT TAKING OR REALLY MISSING!!

ds_med <- ds %>% 
  dplyr::select(SID, VN, fMedsBP) %>% 
  na.omit()

Progression

Generalized Estimating Equations

GEE where predictor is baseline UDBP

# Predictor is baseline UDBP ----------------------------------------

gee_kidney_baseline <- ds %>% 
  prep_mason_data_kidney() %>% 
  mason_geeplot(y = c("lACR", "leGFR"),
            x = "udbpCrBase",
            covars = c("YearsFromBaseline", "ageBase", "Sex", "Ethnicity", 
                       "dmStatus")) %>%
  mason::polish_renaming(rename_gee_kidney)


# GEE table ----------------------------------------------------------

gee_kidney_baseline %>% 
  gee_results_table() %>% 
  pander::pander()

# Plot ---------------------------------------------------------------

plot_gee_results_kidney_base(gee_kidney_baseline,
                 yvars = c("ACR (mg/mmol)", 
                           "eGFR (ml/min/1.73m^2)"))

# Predictor is UDBP at all time points ------------------------------

gee_eGFR_baseline <- ds %>% 
  prep_mason_data_kidney() %>% 
  mason_geeplot(y = c("leGFR"),
            x = "udbpCrBase",
            covars = c("YearsFromBaseline", "ageBase", "Ethnicity", 
                       "DM")) %>%
  mason::polish_renaming(rename_gee_kidney)

# GEE table ----------------------------------------------------------

gee_eGFR_baseline %>% 
  gee_results_table() %>% 
  pander::pander()

# Plot ---------------------------------------------------------------

plot_gee_results_kidney_base(gee_kidney_baseline,
                 yvars = c("ACR (mg/mmol)", 
                           "eGFR (ml/min/1.73m^2)"))

GEE where predictor is UDBP over time

# Predictor is UDBP at all time points ------------------------------

gee_kidney <- ds %>% 
  prep_mason_data_kidney() %>% 
  mason_geeplot(y = c("lACR", "leGFR"),
            x = "udbpCr",
            covars = c("YearsFromBaseline", "ageBase", "Sex", "Ethnicity", 
                       "dmStatus")) %>%
  mason::polish_renaming(rename_gee_kidney)

# GEE table ----------------------------------------------------------

gee_kidney %>% 
  gee_results_table() %>% 
  pander::pander()

# Plot ---------------------------------------------------------------

plot_gee_results_kidney(gee_kidney,
                 yvars = c("ACR (mg/mmol)", 
                           "eGFR (ml/min/1.73m^2)"))

covars <- c("YearsFromBaseline", "ageBase", "Ethnicity", 
                       "DM", "fMedsBP")

gee_terms <- c(
  "<-Xterm",
  "Follow-up Duration (years)",
  "Baseline Age (years)",
  "EthnicityEuropean",
  "Diabetes",
  "fMedsBPYes"
  )

gee_labels <- c(
  "uVDBP:cr (ug/mL)",
  "Follow-up Duration (Years)",
  "Baseline Age (Years)",
  "Ethnicity (European)",
  "Diabetes",
  "Blood Pressure Medication"
  )

# Predictor is UDBP at all time points ------------------------------

gee_eGFR <- ds %>% 
  prep_mason_data_kidney() %>% 
  mason_geeplot(y = c("leGFR"),
            x = "udbpCr",
            covars = covars) %>%
  mason::polish_renaming(rename_gee_kidney)

# GEE table ----------------------------------------------------------

gee_eGFR %>% 
  gee_results_table() %>% 
  pander::pander()

# Plot ---------------------------------------------------------------

plot_gee_results_eGFR(gee_eGFR,
                 yvars = "eGFR (ml/min/1.73m^2)",
                 terms = gee_terms,
                 labels = gee_labels)