R/ddiwas_regression_progress_glm.R

In patrickwu510/ddiwas: What the Package Does (One Line, Title Case)

#' DDIWAS regression with progress bar, regular glm
#'
#' @param regression_table A tibble
#' @param regression_drugs A tibble
#' @return A dataframe \code{ddiwas_results}
#' @export

ddiwas_regression_progress_glm <- function(regression_table, regression_drugs) {

  library(tidyverse)
  library(varhandle)
  library(MASS)

  regression_table <- subset(regression_table, select = -c(person_id))
  regression_drugs<- subset(regression_drugs, select = -c(person_id))

  df.cases <- filter(regression_table, groupc == 1)
  df.controls <- filter(regression_table, groupc == 0)
  nCases <- nrow(df.cases)
  nControls <- nrow(df.controls)

  ddiwas_fx <- function(drug) {

    regression_table["drug"] <- regression_drugs[drug]

    # run model to get statistics
    glmModel <- glm(groupc~drug+obs_length_n+age_n+is_m+is_w, data=regression_table, family=binomial)

    # calculate 2x2 contingency table numbers
    regression_table_cases <- filter(regression_table, groupc == 1)
    regression_table_controls <- filter(regression_table, groupc == 0)
    nA <- sum(regression_table_cases["drug"])
    nB <- sum(regression_table_controls["drug"])
    nC <- nCases-nA
    nD <- nControls-nB

    # store results
    results <- c(drug,
                 glmModel$coefficients[2],
                 sqrt(diag(vcov(glmModel)))[2],
                 glmModel$prob[2],
                 exp(glmModel$coefficients[2]),
                 nA,nB,nC,nD)
  }
  drugs <- names(regression_drugs)
  ddiwas_results <- lapply(drugs, ddiwas_fx)
  ddiwas_results1 <- as.data.frame(do.call(rbind, ddiwas_results))
  names(ddiwas_results1) <- c("rxcui_ingr","coef","se","pval","or","nA","nB","nC","nD")
  ddiwas_results1 <- unfactor(as_tibble(ddiwas_results1))

  return(ddiwas_results1)
}