# 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) 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 <- project_data 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')
# Diet information is not available at baseline. Only VN 3 data is available at the present time. tab_nums("subchar_vitd_baseline", "Subject characteristics across vitamin D status at baseline") subchar_vitd_base <- tableone::CreateTableOne( vars = c("Age", "Sex", "Ethnicity", "BMI", "Waist", "eGFR", "ACR", "UrineCreatinine", "UrineMicroalbumin", "UrinaryCalcium", "UDBP", "udbpCr", "VitaminD", "PTH", "MET", "Systolic", "Diastolic", "MeanArtPressure", "diet_milk", "diet_cal", "diet_supp_cal", "diet_supp_vitd", "dmStatus"), strata = c("vitdStatus"), data = dsBase, factorVars = c("Sex", "Ethnicity", "dmStatus") ) %>% print(nonnormal = c("UDBP", "ACR", "UrineMicroalbumin"), quote = FALSE, noSpaces = TRUE) %>% knitr::kable()
# Diet information is not available at baseline. Only VN 3 data is available at the present time. tab_nums("subchar_vitd", "Subject characteristics across visit numbers") # for sensitivity analysis m <- ds %>% dplyr::filter(!is.na(VitaminD)) %>% dplyr::filter(!VN == "6") subchar_vitd_base <- tableone::CreateTableOne( vars = c("Age", "Sex", "Ethnicity", "BMI", "Waist", "MET", "VitaminD", "PTH", "ACR", "eGFR", "UDBP", "udbpCr", "Systolic", "Diastolic", "MeanArtPressure", "dmStatus"), strata = c("VN"), data = m, factorVars = c("Sex", "Ethnicity", "dmStatus") ) %>% print(nonnormal = c("UDBP", "udbpCr", "ACR"), quote = FALSE, noSpaces = TRUE) %>% knitr::kable() ds %>% dplyr::select(VN, PTH) %>% dplyr::group_by(VN) %>% na.omit() %>% summarise(n = n())
# Clean data vitd <- dsBase %>% dplyr::mutate(vitdLabels = ifelse(vitdStatus == "Sufficient", "Sufficient \n(≥75 nmol/L)", ifelse(vitdStatus == "Insufficient", "Insufficient \n(50 - 75 nmol/L)", "Deficient \n(<50 nmol/L)")), vitdLabels = factor(vitdLabels, levels = c("Sufficient \n(≥75 nmol/L)", "Insufficient \n(50 - 75 nmol/L)", "Deficient \n(<50 nmol/L)"), ordered = TRUE)) %>% select(vitdStatus, vitdLabels, udbpCr, UDBP) %>% na.omit() # Box plot of uVDBP in different albuminuria categories vitd %>% box_plot_slides("vitdLabels", "log(udbpCr)", "Vitamin D Status", "log uVDBP:cr (μg/mmol)") # value vitd %>% dplyr::group_by(vitdStatus) %>% dplyr::summarise(n = n(), Median = median(log(udbpCr)), lower = quantile(log(udbpCr), probs = 0.25), higher = quantile(log(udbpCr), probs = 0.75), IQR = IQR(log(udbpCr)) ) # n vitd %>% group_by(vitdStatus) %>% summarise(n = n()) # ANOVA anova <- aov(formula = log(udbpCr)~vitdStatus, data = vitd) summary(anova) TukeyHSD(anova) rm(anova)
# Scatterplot of ACR and uVDBP ---------------------------------- dsBase %>% # dplyr::filter(acrStatus != "Normoalbuminuria") %>% scatter_plot_transparent("log(udbpCr)", "log(VitaminD)", "log uVDBP:cr (μg/mmol)", "log serum 25(OH)D (nmol/L)") # Spearman Correlation ------------------------------------------ dsBase %>% cor.test(formula = ~ log(udbpCr) + log(VitaminD), data = ., method = "spearman") # Linear Regression --------------------------------------------- dsBase %>% dplyr::mutate(lVitD = log(VitaminD)) %>% prep_mason_data() %>% mason_glm(y = "logudbpCr", x = "lVitD", covars = c("ageBase", "Sex", "Ethnicity", "MET", "Season", "dmStatus") ) %>% gee_results_table() %>% pander::pander()
# Scatterplot dsBase %>% scatter_plot("PTH", "VitaminD", "Parathyroid Hormone (pmol/L)", "Serum 25(OH)D (nmol/L)") # Linear Regression --------------------------------------------- dsBase %>% prep_mason_data() %>% mason_glm(y = "VitaminD", x = "PTH" ) %>% dplyr::filter(!term == "(Intercept)") %>% dplyr::mutate(p = round(p.value, 2), p = ifelse(p == "0", "<0.001", p), estCI = paste0(round(estimate, 2), " (", round(conf.low, 2), ", ", round(conf.high, 2), ")")) %>% dplyr::select(Yterms, Xterms, term, estCI, p) %>% # tidyr::spread(Yterms, estCI) %>% pander::pander()
# Box plot of 25(OH)D in different seasons dsBase %>% box_plot("Season", "VitaminD", "Seasons", "Serum 25(OH)D (nmol/L)") # n dsBase %>% group_by(Season) %>% summarise(n = n()) # ANOVA anova <- aov(formula = VitaminD~Season, data = dsBase) summary(anova) TukeyHSD(anova) rm(anova)
ds %>% plot_progress(yvar = "VitaminD", ylab = "Serum 25(OH)D") # Complete data dsComplete %>% plot_progress(yvar = "VitaminD", ylab = "Serum 25(OH)D") # ANOVA ------------------------------------------------------------ anova <- aov(formula = VitaminD~fVN, data = dsComplete) summary(anova) TukeyHSD(anova) rm(anova)
# Predictor is baseline UDBP ---------------------------------------- gee_vitd_baseline <- ds %>% prep_mason_data_vitd() %>% mason_geeplot(y = "lVitD", x = "udbpCrBase", covars = c("YearsFromBaseline", "ageBase", "Sex", "Ethnicity", "BMI", "MET", "dmStatus", "Season"), intvar = "BMI") %>% mason::polish_renaming(rename_gee_vitd) # GEE table ---------------------------------------------------------- gee_vitd_baseline %>% gee_results_table() %>% # tidyr::spread(Yterms, estCI) %>% pander::pander() # Plot --------------------------------------------------------------- plot_gee_results_vitd_base(gee_vitd_baseline, yvars = "Serum 25(OH)D (nmol/L)") # Plot interaction --------------------------------------------------- # dsBase %>% # interac_plot(yvar = "PTH", # intvar = "Season", # ylab = "Parathyroid Hormone (pmol/L)") #Serum 25(OH)D (nmol/L) #Parathyroid Hormone (pmol/L) # # rms::ols(VitaminD ~ udbpCr * Season, data = dsBase)
# Predictor is UDBP ---------------------------------------- gee_vitd <- ds %>% prep_mason_data_vitd() %>% mason_geeplot(y = "lVitD", x = "udbpCr", covars = c("YearsFromBaseline", "ageBase", "Sex", "Ethnicity", "BMI", "MET", "dmStatus", "Season"), intvar = "BMI") %>% mason::polish_renaming(rename_gee_vitd) # GEE table ---------------------------------------------------------- gee_vitd %>% gee_results_table() %>% # tidyr::spread(Yterms, estCI) %>% pander::pander() # Plot --------------------------------------------------------------- plot_gee_results_vitd(gee_vitd, yvars = "Serum 25(OH)D (nmol/L)") # Plot interaction --------------------------------------------------- # ds %>% # interac_plot(yvar = "VitaminD", # intvar = "Season", # ylab = "Serum 25(OH)D (nmol/L)") # # rms::ols(VitaminD ~ udbpCr * Season, data = ds)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.