REMA-4 Données manquantes

Jul 1, 2020 12 min read

Précédemment nous avons vu différents moyens d’importer, manipuler, visualiser et analyser nos données. Je vous ai souvent mentionné que les données manquantes étaient parfois problématiques.

Maintenant on n’a pas le choix, il faut regarder cela et y voir de plus près…

Avoir un modèle qui roule c’est bien mais on doit savoir quels sont les limites et difficultés à l’interpréter en présence de données manquantes.

Encore une fois, par souci de transparence il est souhaitable d’indiquer ce que l’on fait avec les données manquantes lorsqu’on en a dans notre BD.

Parfois on peut bien chercher mais cela ne sert à rien, les données manquantes sont là pour rester…

source

Bien sûr, je ne réinvente rien! Je vous réfère donc au livre de S VanBuuren : Flexible imputation for missing data. Je me suis basé sur la lecture de ce livre pour avancer moi-même dans ma compréhension des données manquantes et donc autant vous faire profiter de mes lectures et réflexions.

Parmi les concepts, il n’existe pas de seuil de magique pour dire à partir de quel pourcentage de NA on s’inquiète. Il faut plus connaitre l’impact des NA et leur origine. Cela va varier beaucoup selon la BD et l’analyse/patron d’absence des données.

Commençons maintenant par charger certains packages nécessaires au tutoriel aujourd’hui.

library(tidyverse)

## -- Attaching packages -------------------------------- tidyverse 1.3.0 --

## v ggplot2 3.3.2     v purrr   0.3.4
## v tibble  3.0.1     v dplyr   1.0.0
## v tidyr   1.1.0     v stringr 1.4.0
## v readr   1.3.1     v forcats 0.5.0

## -- Conflicts ----------------------------------- tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()

library(readxl)
library(funModeling)

## Le chargement a nécessité le package : Hmisc

## Le chargement a nécessité le package : lattice

## Le chargement a nécessité le package : survival

## Le chargement a nécessité le package : Formula

## 
## Attachement du package : 'Hmisc'

## The following objects are masked from 'package:dplyr':
## 
##     src, summarize

## The following objects are masked from 'package:base':
## 
##     format.pval, units

## funModeling v.1.9.4 :)
## Examples and tutorials at livebook.datascienceheroes.com
##  / Now in Spanish: librovivodecienciadedatos.ai

1. Généralités sur les données manquantes

Dans R, le format standard pour les spécifier est la notation NA. C’est le seul format officiel.

Il va donc être bien important de s’assurer que le codage de ces données répond à ce format lorsque nous allons travailler avec ces dernières.

1.1 Spécifier que les données manquantes sont sous un format particulier

Prenons un exemple concrêt. Je simule un vecteur de 10 données dont 3 sont manquantes

sim_na0 <- c(1,  2, na, NA, 2, 3, 3, 2, 1, Na) #ici on doit noter que si je ne mets pas de parenthèses, mis à part le format NA les "na" et "Na" ne sont pas reconnus donc j'au un message d'erreur

## Error in eval(expr, envir, enclos): objet 'na' introuvable

#je dois donc spécifier que ce sont des characters en mettant des guillemets.
sim_na <- c(1,  2, "na", NA, 2, 3, 3, 2, 1, "Na")
class(sim_na)

## [1] "character"

J’obtiens alors un vecteur qui par défaut va être reconnu comme étant formé de character. Pour changer tout cela il faut que je change dans un premier temps la façon dont sont entrées ces données manquantes afin qu’elles soient reconnues comme NA. Pour cela rien de plus simple: je vais utiliser la fonction na_if du paquet dplyr.

NB: je vous donne uniquement une méthode, il en existe plusieurs et peut être des plus efficaces…

sim_na %>% na_if("na") %>% na_if("Na")

##  [1] "1" "2" NA  NA  "2" "3" "3" "2" "1" NA

Donc on a bien corrigé le tout et changé notre vecteur en faisant reconnaitre nos NA comme telles.

Cependant, comme vous le voyez ci-dessous par défaut R a donc codé le sim_na comme character du fait que les formats différents de données manquantes étaient de ce format.

class(sim_na)

## [1] "character"

On doit donc retransformer le vecteur de cette façon:

sim_na=as.double(sim_na)

## Warning: NAs introduits lors de la conversion automatique

sim_na

##  [1]  1  2 NA NA  2  3  3  2  1 NA

COOL !!!…

1.2 Connaitre le nombre de données manquantes d’une variable

Comme précédemment vu dans les différents tutoriels, la fonction summary() fait un excellent travail pour donner la base lorsque la variable est numérique

summary(sim_na)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
##     1.0     1.5     2.0     2.0     2.5     3.0       3

Si nous avons une variable non numérique

sim2 = c("Rouge", "Blanc", "Bleu", NA, "Vert", "Jaune")
summary(sim2)

##    Length     Class      Mode 
##         6 character character

Dans ce cas là la fonction summary ne marche pas. On utilise quelque chose de Beaucoup plus simple. Nous voulons savoir combien de fois la réponse is.na() est vraie !!!

is.na(sim2) # applique la question à toute mes lignes

## [1] FALSE FALSE FALSE  TRUE FALSE FALSE

sum(is.na(sim2)) # compte le tout...

## [1] 1

Notez bien que la dernière façon sum(is.na()) fonctionne quelle que soit le type de variable…

2. Comment gérer les données manquantes

2.1 Omission

Simplement!!!

Source

Je sais c’est fatigant mais c’est ce que font par défaut les logiciels classiques. Cela va mener à une analyse de cas complets (complete case analysis) mais souvent c’est fait de façon silencieuse sans vous le mentionner dans la plupart des logiciels. Avec R, c’est différent (et parfois crispant) car lui va tout le temps vous le mettre sous le nez si vous avez des NA. Vous allez devoir expliciter les choix que vous prenez pour tenir compte ou pas de ces dernières.

2.1.1 Ne pas tenir compte des données manquantes

Un exemple concret maintenant.

Imaginez que vous avez collecté des données d’intervalle velâge-vêlage dans un troupeau sur 6 vaches (Un petit troupeau).

IVV <- c(355, 388, 400, NA, 377, 380, 480, 420)
mean(IVV)

## [1] NA

Donc R ne peut pas me répondre. Il considère que comme toutes les données ne sont pas complètes la réponse n’a pas de sens…

Je dois indiquer que je veux ne pas tenir compte des données manquantes.

mean(IVV, na.rm=TRUE) # j'active la commande remove NA

## [1] 400

2.1.2 Application dans une base de données

Comme précédemment je vais utiliser la même base de données concernant les veaux du livre Veterinary Epidemiology Research de Ian Dohoo.

veaux <- read_excel("calf.xlsx") 
head(veaux, n=4)

## # A tibble: 4 x 14
##    case   age breed   sex  attd  dehy   eye  jnts  post pulse  resp  temp   umb
##   <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1  1670     5     1     1     2  12      NA     0     2    NA    NA  37.6     0
## 2  8124     3     2     0     1  13.5     0     0     0   130   120  39.2     0
## 3  6954     2     3     1     2  NA       1     0     2    NA    NA  NA       1
## 4  2737     3     4     1     1   5       0     0     0   132    40  38.6     0
## # ... with 1 more variable: sepsis <dbl>

#tout d'abord je veux prédire le sepsis comme précédemment

class(veaux$sepsis) # d'emblée elle est reconnue comme numérique je dois donc m'assurer que c'est un facteur

## [1] "numeric"

veaux <- veaux %>% mutate(septic = as_factor(sepsis)) # variable à prédire

#ensuite je crée la variable qui sépare mon age en quintiles
veaux <- veaux %>% mutate(age_c = factor(ntile(age, 5)))

cross_plot(veaux, input = "age_c", target = "septic", plot_type="percentual")

Pour l’instant faisons fi des messages de warning.

Regardons les données manquantes de cette base de données.

colSums(is.na(veaux))

##   case    age  breed    sex   attd   dehy    eye   jnts   post  pulse   resp 
##      0      1      0      2      6     14     18     11      6      9     19 
##   temp    umb sepsis septic  age_c 
##      7     11      0      0      1

Ici je vois donc que j’ai différentes données manquantes dans différentes colonnes.

Je pourrais aussi le visualiser plus facilement en utilisant un package plus spécifiquement utile aux données manquantes.

C’est le package naniar.

library(naniar)

vis_miss(veaux, sort_miss = TRUE) #je classifie aussi par le nombre de données manquantes

On peut voir ici que j’ai 2.6% de données manquantes (moyenne sur l’ensemble de ma base de données). La composante respiratoire est la plus fréquemment affectée avec 7.48% de données qui manquent…

Je peux aussi regarder si des patrons de données manquantes sont observés plus fréquemment que d’autres avec la fonction gg_miss_upset qui est une fonction qui me permet de regarder les combinaisons de données manquantes les plus fréquentes.

gg_miss_upset(veaux)

Grâce à cette figure je vois le type d’associations de patrons de données manquantes que j’ai dans mon jeu de données ainsi que le nombre de ces associations de données manquantes. Par exemple, j’ai 3 cas où je trouve à la fois des NA chez le même animal pour les variables jnts, umb, dehy, eye et resp.

Je peux aussi tout simplement visualiser les données manquantes par variable classées par ordre décroissant sous la forme d’une figure lollypop avec la fonction gg_miss_var().

gg_miss_var(veaux)

Un autre package qui est très utilisé pour les données manquantes est le package Amelia

library(Amelia)

## Le chargement a nécessité le package : Rcpp

## ## 
## ## Amelia II: Multiple Imputation
## ## (Version 1.7.6, built: 2019-11-24)
## ## Copyright (C) 2005-2020 James Honaker, Gary King and Matthew Blackwell
## ## Refer to http://gking.harvard.edu/amelia/ for more information
## ##

missmap(veaux,
 main = "Autre figure avec le paquet amelia",
 col=c("red", "skyblue"),
 y.labels = NULL,
 y.at = NULL)

## Warning: Unknown or uninitialised column: `arguments`.

## Warning: Unknown or uninitialised column: `arguments`.

Notez que dans ce tutoriel je ne parle que de visualisation des données manquantes de votre BD. Je ne parle pas des possibilités de gérer ces données manquantes (ce que vous permettrons ces mêmes paquets).

Si vous avez peu de données manquantes (peu?) alors on peut simplement demander à R de les ignorer par la commande na.action où vous spécifiez ce que vous voulez faire en disant que vous voulez omettre les données manquantes.

model <- glm(septic~age_c + temp, family = binomial(link = "logit"), data=veaux, na.action=na.omit)
summary(model)

## 
## Call:
## glm(formula = septic ~ age_c + temp, family = binomial(link = "logit"), 
##     data = veaux, na.action = na.omit)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.3644  -0.8549  -0.6401   1.2137   2.1423  
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)   
## (Intercept)  8.10264    3.30192   2.454  0.01413 * 
## age_c2      -0.26252    0.41940  -0.626  0.53135   
## age_c3      -1.39388    0.48008  -2.903  0.00369 **
## age_c4      -0.91883    0.44598  -2.060  0.03938 * 
## age_c5      -1.08348    0.46295  -2.340  0.01927 * 
## temp        -0.21923    0.08595  -2.551  0.01075 * 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 295.65  on 245  degrees of freedom
## Residual deviance: 278.67  on 240  degrees of freedom
##   (8 observations deleted due to missingness)
## AIC: 290.67
## 
## Number of Fisher Scoring iterations: 4

Je pourrai même spécifier à R de le faire par défaut pour tout mon code en précisant l’option.

Notez que contrairement aux autres logiciels, R ne cache pas les NA par défaut. Il faut que nous le fassions manuellement. Parfois c’est C###!!!! mais cela permet de bien comprendre la structure de sa BD toujours à des fins de reproductibilité.

options(na.action=na.omit)

Mais admettons maintenant que je veuille lier mes prédictions à ma base de données (pour voir par exemple le lien entre les prédictions et les observations).

#j'utilise la fonction cbind qui me permet de lier 2 BD  
prediction <- cbind(veaux, predict(model))

## Error in data.frame(..., check.names = FALSE): les arguments impliquent des nombres de lignes différents : 254, 246

Ici cela ne fonctionne pas car les dimensions de mes 2 BD sont différentes concernant le nombre d’observations (puisque par défaut mon modèle logistique ne contient pas les données manquantes).

naprint(na.action(model))

## [1] "8 observations deleted due to missingness"

Je dois donc préciser à R que je veux qu’il omette les NA de ma BD (pour les variables incluses dans le modèle afin de pouvoir ensuite les joindre au modèle (par l’argument cbind qui joint des colonnes).

veaux2 <- cbind(na.omit(veaux[, c("septic", "age_c", "temp")]),
                     predicted = predict(model))

Ensuite je peux donc travailler avec l’ensemble des données de mon modèle (celles qui ont servies à son élaboration).

Mais cela dépasse le cadre de ce tutoriel

2.2 Imputation des données manquantes

La principale limite à l’élimination de données manquantes est donc qu’on perd énormément d’information surtout si les donnnées manquantes sont réparties dans toute la base de données.

Prenons un exemple concrêt: Imaginons qu’on ait 100 observations (1 observation= 1 ligne avec \(Y\) et 5 covariables \(x\)) et 20 données manquantes parmi les différents prédicteurs \(x\). Si certaines de ces données manquantes concernent la même ligne je pourrais avoir un nombre de lignes incomplètes <20 (par exemple dans le meilleur des cas (si on est TRÈS chanceux) je pourrai avoir uniquement 4 individus pour lesquels je n’ai aucun des 5 prédicteurs \(x\) (j’ai juste leur issue observée \(Y\)). Ces 4 individus correspondraient à mes 20 données manquantes (4*5)

À ce moment là cela me fait uniquement 4% de perte si j’ignore mes données manquantes. Par contre si mes données manquantes sont réparties dans des lignes différentes à ce moment là je perds 20% de mes observations pour une analyse en multivarié avec des cas complets. Ici cela peut poser de sérieux problèmes surtout concernant la possible relation/équivalence entre ce modèle et la réalité qu’on vise à estimer, surtout si le fait d’avoir ces NA dépend d’une autre particularité de nos données…

Imaginez que tous les veaux septiques sont ceux où on n’a pas d’info sur un paramètre parce qu’il est trop altéré ou difficile à évaluer. A ce moment votre modèle risque d’être impacté par ces spécificités et de ne pas mettre en lumière ce paramètre alors qu’il est très important.

Je n’élaborerai pas trop là dessus mais c’est donc TRÈS important de se demander quel est le patron de ces données qui manquent et aussi de savoir comment on fait pour imputer ces données soit par des méthodes simples (ex la moyenne ou la prédiction d’une imputation). Pour cela je vous réfère aux documents cités précédemment.

Fin de cette introduction

Pour toute information ou problème me contacter s.buczinski@umontreal.ca

Sébastien Buczinski

Professeur titulaire clinique ambulatoire bovine

Ma recherche se focalise sur la santé des veaux, les maladies respiratoires bovines, les tests diagnostiques en l’absence de test de référence parfaits et la médecine vétérinaire factuelle.