R/polycumul.R

In INESSSQC/polymedic: Calculs de la polypharmacie

#' Polypharmacie Cumulée
#'
#' La polypharmacie cumulée considère qu'il y a eu consommation d'un médicament s'il y a au moins une journée dans l'intervalle [t1, t2] un médicament pris par la personne.
#'
#' \strong{\code{period}} : Divise l'intervalle d'analyse en \code{period} parties égales. Toutefois, le nombre de jours est arrondi. Par exemple, si un intervalle d'analyse de 10 jours est diviser en trois parties, il y aura une section de 3 jours -> [0,3[, de 4 jours -> [3, 7[ et de 3 jours -> [7, 10[.
#'
#' @param x Base de données créée par \code{polyBD()}.
#' @param t1,t2 "AAAA-MM-JJ". Intervalle de temps [t1, t2] servant d'analyse.
#' @param period Par défaut 4. INTEGER. Divise l'intervalle d'analyse [\code{t1}, \code{t2}] en \code{period} partie(s) égale(s) (Voir détails). Le résultat de la consommation du nombre de médicament est la moyenne des périodes.
#' @param catage Borne(s) inférieure(s) des groupes d'âge.
#'
#' @import dplyr
#' @import data.table
#'
#' @export
polycumul <- function(x, t1, t2, period = 4, catage = c(60,70,80,90)){

  ## Noms des colonnes
  # Conserver
  xcols_ <- names(x)
  # Modifier
  names(x)[[1]] <- "BenBanls"
  names(x)[[2]] <- "DIN"
  # Liste des ID
  xID <- unique(x$BenBanls)

  # Convertir t1 et t2 en numeric
  t1 <- as.Date(t1)
  t2 <- as.Date(t2)
  # Intervalles des périodes
  bornes <- lapply(0:period, function(x) round(t1 + ((t2-t1) / period) * x))

  # Consommation par période
  cumul <- data.table()  # tableau à insérer résultats
  for(i in seq(period)){
    # Médicaments consommé dans la période
    dt <- x[PerDebut <= bornes[[i+1]] & PerFinAjust >= bornes[[i]]]
    # 1 observation par BenBanls + DIN
    dt <- unique(dt[, .(BenBanls, DIN)])
    # Nombre de médicaments (DIN) par BenBanls
    dt <- dt[, .(n = .N), .(BenBanls)]

    # Ajouter les BenBanls manquants - sans consommation
    xID_dt <- xID[!xID %in% dt$BenBanls]
    # Tableau indiquant leur consommation = 0
    if(length(xID_dt) != 0){
      xID_dt <- data.table(BenBanls = xID_dt,
                           n = 0L)
      # Ajout des observations de xID_dt à dt
      dt <- rbind(dt, xID_dt)
    }
    # Afficher période
    dt[, period := i]
    # Ajout des observations de dt à cumul
    cumul <- rbind(cumul, dt)
  }
  # Tri cumul par personne
  setorder(cumul, BenBanls)
  # Moyenne de la consommation

  ## Nombre d'usager
  nID <- length(xID)


  ## Résultats
  # BD d'analyse
  cumul_BD <- cumul[, Mean := mean(n), .(BenBanls)]
  # Statistiques
  cumul_stats <- unique(cumul_BD[, .(BenBanls, Mean)])

  cumul_stats_pop <- cumul_stats[, .(Moyenne = mean(Mean),
                                     `Écart-type` = sd(Mean),
                                     Min = min(Mean),
                                     P5 = quantile(Mean, 0.05),
                                     P10 = quantile(Mean, 0.10),
                                     Q1 = quantile(Mean, 0.25),
                                     `Médiane` = median(Mean),
                                     Q3 = quantile(Mean, 0.75),
                                     P90 = quantile(Mean, 0.90),
                                     P95 = quantile(Mean, 0.95),
                                     Max = max(Mean),
                                     n = nID)]

  # Ajout de la date de naissance
  cumul_stats <- merge(cumul_stats, unique(attributes(x)$Naissance),
                       by.x = "BenBanls", by.y = xcols_[[1]])
  # Convertir date de naissance en Age (selon la date d'analyse de départ)
  cumul_stats[
    , Age := as.integer(floor(difftime(t1, DatNais, units = "days") / 365.25))
  ]
  # Inclure tous les âges dans 'catage'
  catage <- c(catage, max(cumul_stats[, Age]) + 1)
  # Création de la catégorie d'age
  cumul_stats[
    , `:=` (CatAge = cut(Age, catage, right = FALSE),
            # Supprimer les colonnes non nécessaires
            Age = NULL,
            DatNais = NULL)
  ]
  # Ajout du Sexe
  cumul_stats <- merge(cumul_stats, unique(attributes(x)$Sexe),
                       by.x = "BenBanls", by.y = xcols_[[1]])
  setnames(cumul_stats, "BenBanls", xcols_[[1]])

  # Statistiques descriptives
  # Par Age et Sexe
  nAS <- cumul_stats[, .(n = .N), keyby = .(CatAge, Sexe)]  # Nombre d'observation par groupe
  cumul_stats_AS <- cumul_stats[, .(Moyenne = mean(Mean),
                                    `Écart-type` = sd(Mean),
                                    Min = min(Mean),
                                    P5 = quantile(Mean, 0.05),
                                    P10 = quantile(Mean, 0.10),
                                    Q1 = quantile(Mean, 0.25),
                                    `Médiane` = median(Mean),
                                    Q3 = quantile(Mean, 0.75),
                                    P90 = quantile(Mean, 0.90),
                                    P95 = quantile(Mean, 0.95),
                                    Max = max(Mean)),
                                keyby = .(CatAge, Sexe)]
  cumul_stats_AS <- merge(cumul_stats_AS, nAS)
  # Par Age
  nA <- cumul_stats[, .(n = .N), keyby = .(CatAge)]  # Nombre d'observation par groupe
  cumul_stats_A <- cumul_stats[, .(Moyenne = mean(Mean),
                                    `Écart-type` = sd(Mean),
                                    Min = min(Mean),
                                    P5 = quantile(Mean, 0.05),
                                    P10 = quantile(Mean, 0.10),
                                    Q1 = quantile(Mean, 0.25),
                                    `Médiane` = median(Mean),
                                    Q3 = quantile(Mean, 0.75),
                                    P90 = quantile(Mean, 0.90),
                                    P95 = quantile(Mean, 0.95),
                                    Max = max(Mean)),
                                keyby = .(CatAge)]
  cumul_stats_A <- merge(cumul_stats_A, nA)
  # Par Sexe
  nS <- cumul_stats[, .(n = .N), keyby = .(Sexe)]  # Nombre d'observation par groupe
  cumul_stats_S <- cumul_stats[, .(Moyenne = mean(Mean),
                                   `Écart-type` = sd(Mean),
                                   Min = min(Mean),
                                   P5 = quantile(Mean, 0.05),
                                   P10 = quantile(Mean, 0.10),
                                   Q1 = quantile(Mean, 0.25),
                                   `Médiane` = median(Mean),
                                   Q3 = quantile(Mean, 0.75),
                                   P90 = quantile(Mean, 0.90),
                                   P95 = quantile(Mean, 0.95),
                                   Max = max(Mean)),
                               keyby = .(Sexe)]
  cumul_stats_S <- merge(cumul_stats_S, nS)



  # Distribution des maximums
  dist_max <- unique(cumul_BD[, .(BenBanls, Mean)])
  dist_max[, Mean := round(Mean)]
  seq_values <- 0:max(dist_max[, Mean]) #toutes les valeurs à afficher -> 0 à valeur max
  dist_max <- dist_max[, .(Freq = .N), Mean] #fréquence pour chaque valeur des maximums
  names(dist_max) <- c("Valeur", "Freq")
  if(length(seq_values[!seq_values %in% dist_max$Valeur]) != 0){
    dist_max <- rbind(dist_max, data.table(Valeur = seq_values[!seq_values %in% dist_max$Valeur],
                                           Freq = 0)) #ajouter les valeurs manquante -> Freq = 0
  }
  setorder(dist_max, Valeur) #tri croissant par Max
  dist_max[, Pourcentage := Freq / sum(Freq) * 100]  # pourcentage
  dist_max[, Cumul := cumsum(Pourcentage)] #pourcentage cumulé


  # Liste des résultats
  cumul <- list(BD = cumul_BD,  # BD d'analyse
                stats = cumul_stats_pop,  # Statistiques
                stats_group = list(`CatAge & Sexe` = cumul_stats_AS,
                                   CatAge = cumul_stats_A,
                                   Sexe = cumul_stats_S),
                dist_max = dist_max)  # Distribution des maximums

}