1 Introduction

A two-step approach to imputing missing data in metabolomics. Step 1 uses a random forest classifier to classify missing values as either Missing Completely at Random/Missing At Random (MCAR/MAR) or Missing Not At Random (MNAR). MCAR/MAR are combined because it is often difficult to distinguish these two missing types in metabolomics data. Step 2 imputes the missing values based on the classified missing mechanisms, using the appropriate imputation algorithms. Imputation algorithms tested and available for MCAR/MAR include Bayesian Principal Component Analysis (BPCA), Multiple Imputation No-Skip K-Nearest Neighbors (Multi_nsKNN), and Random Forest. Imputation algorithms tested and available for MNAR include nsKNN and a single imputation approach for imputation of metabolites where left-censoring is present.

2 Installation

The following code chunk depicts how to install MAI from Bioconductor

if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")

BiocManager::install("MAI")

3 Using MAI when your data is a data.frame or matrix

# Load the MAI package
library(MAI)
# Load the example data with missing values
data("untargeted_LCMS_data")
# Set a seed for reproducibility
## Estimating pattern of missingness involves imposing MCAR/MAR into the data
## these are done at random and as such may slightly change the results of the
## estimated parameters.
set.seed(137690)
# Impute the data using BPCA for predicted MCAR value imputation and
# use Single imputation for predicted MNAR value imputation

Results = MAI(data_miss = untargeted_LCMS_data, # The data with missing values
          MCAR_algorithm = "BPCA", # The MCAR algorithm to use
          MNAR_algorithm = "Single", # The MNAR algorithm to use
          assay_ix = 1, # If SE, designates the assay to impute
          forest_list_args = list( # random forest arguments for training
            ntree = 300,
            proximity = FALSE
        ),
          verbose = TRUE # allows console message output
        )
## Estimating pattern of missingness
## Imposing missingness
## Generating features
## Training
## Predicting
## Imputing
# Get MAI imputations
Results[["Imputed_data"]][1:5, 1:5] # show only 5x5
##          [,1]     [,2]     [,3]     [,4]     [,5]
## [1,] 12.87164 12.70516 12.11793 12.07897 11.18725
## [2,] 12.10813 12.36043 12.08463 12.07897 14.30504
## [3,] 12.10813 11.16815 12.08463 12.07897 11.79195
## [4,] 12.01395 12.28791 12.12017 12.07897 10.05397
## [5,] 11.32825 12.55314 12.23772 11.31292 12.28869
# Get the estimated mixed missingness parameters
Results[["Estimated_Params"]]
## $Alpha
## [1] 20
## 
## $Beta
## [1] 80
## 
## $Gamma
## [1] 60

These parameters estimate the ratio of MCAR/MAR to MNAR in the data. The parameters \(\alpha\) and \(\beta\) separate high, medium, and low average abundance metabolites, while the parameter \(\gamma\) is used to impose missingness in the medium and low abundance metabolites. A smaller \(\alpha\) corresponds to more MCAR/MAR being present, while larger \(\beta\) and \(\gamma\) values imply more MNAR values being present. The returned estimated parameters are then used to impose known missingness in the complete subset of the input data. Subsequently, a random forest classifier is trained to classify the known missingness in the complete subset of the input data. Once the classifier is established it is applied to the unknown missingness of the full input data to predict the missingness. Finally, the missing values are imputed using a specific algorithm, chosen by the user, according to the predicted missingness mechanism.

4 Using MAI when your data is a SummarizedExperiment (SE) class

# Load the SummarizedExperiment package
suppressMessages(
  library(SummarizedExperiment)
)

# Load the example data with missing values
data("untargeted_LCMS_data")
# Turn the data to a SummarizedExperiment
se = SummarizedExperiment(untargeted_LCMS_data)
# Set a seed for reproducibility
## Estimating pattern of missingness involves imposing MCAR/MAR into the data
## these are done at random and as such may slightly change the results of the
## estimated parameters.
set.seed(137690)
# Impute the data using BPCA for predicted MCAR value imputation and
# use Single imputation for predicted MNAR value imputation

Results = MAI(se, # The data with missing values
          MCAR_algorithm = "BPCA", # The MCAR algorithm to use
          MNAR_algorithm= "Single", # The MNAR algorithm to use
          assay_ix = 1, # If SE, designates the assay to impute
          forest_list_args = list( # random forest arguments for training
            ntree = 300,
            proximity = FALSE
        ),
          verbose = TRUE # allows console message output
        )
## Estimating pattern of missingness
## Imposing missingness
## Generating features
## Training
## Predicting
## Imputing
# Get MAI imputations
assay(Results)[1:5, 1:5] # show only 5x5
##          [,1]     [,2]     [,3]     [,4]     [,5]
## [1,] 12.87164 12.70516 12.11793 12.07897 11.18725
## [2,] 12.10813 12.36043 12.08463 12.07897 14.30504
## [3,] 12.10813 11.16815 12.08463 12.07897 11.79195
## [4,] 12.01395 12.28791 12.12017 12.07897 10.05397
## [5,] 11.32825 12.55314 12.23772 11.31292 12.28869
# Get the estimated mixed missingness parameters
metadata(Results)[["meta_assay_1"]][["Estimated_Params"]]
## $Alpha
## [1] 20
## 
## $Beta
## [1] 80
## 
## $Gamma
## [1] 60

5 Session Information

sessionInfo()
## R Under development (unstable) (2023-10-22 r85388)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 22.04.3 LTS
## 
## Matrix products: default
## BLAS:   /home/biocbuild/bbs-3.19-bioc/R/lib/libRblas.so 
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.0
## 
## locale:
##  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
##  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
##  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
##  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
## 
## time zone: America/New_York
## tzcode source: system (glibc)
## 
## attached base packages:
## [1] stats4    stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] SummarizedExperiment_1.33.0 Biobase_2.63.0             
##  [3] GenomicRanges_1.55.0        GenomeInfoDb_1.39.0        
##  [5] IRanges_2.37.0              S4Vectors_0.41.0           
##  [7] BiocGenerics_0.49.0         MatrixGenerics_1.15.0      
##  [9] matrixStats_1.0.0           caret_6.0-94               
## [11] lattice_0.22-5              ggplot2_3.4.4              
## [13] MAI_1.9.0                   BiocStyle_2.31.0           
## 
## loaded via a namespace (and not attached):
##  [1] bitops_1.0-7            pROC_1.18.4             rlang_1.1.1            
##  [4] magrittr_2.0.3          e1071_1.7-13            compiler_4.4.0         
##  [7] vctrs_0.6.4             reshape2_1.4.4          stringr_1.5.0          
## [10] pkgconfig_2.0.3         crayon_1.5.2            fastmap_1.1.1          
## [13] XVector_0.43.0          utf8_1.2.4              rmarkdown_2.25         
## [16] prodlim_2023.08.28      itertools_0.1-3         purrr_1.0.2            
## [19] xfun_0.40               randomForest_4.7-1.1    zlibbioc_1.49.0        
## [22] cachem_1.0.8            jsonlite_1.8.7          recipes_1.0.8          
## [25] DelayedArray_0.29.0     parallel_4.4.0          R6_2.5.1               
## [28] bslib_0.5.1             stringi_1.7.12          parallelly_1.36.0      
## [31] rpart_4.1.21            lubridate_1.9.3         jquerylib_0.1.4        
## [34] Rcpp_1.0.11             bookdown_0.36           iterators_1.0.14       
## [37] knitr_1.44              future.apply_1.11.0     Matrix_1.6-1.1         
## [40] splines_4.4.0           nnet_7.3-19             timechange_0.2.0       
## [43] tidyselect_1.2.0        abind_1.4-5             yaml_2.3.7             
## [46] timeDate_4022.108       doParallel_1.0.17       codetools_0.2-19       
## [49] doRNG_1.8.6             listenv_0.9.0           tibble_3.2.1           
## [52] plyr_1.8.9              withr_2.5.1             evaluate_0.22          
## [55] future_1.33.0           survival_3.5-7          proxy_0.4-27           
## [58] pillar_1.9.0            BiocManager_1.30.22     rngtools_1.5.2         
## [61] foreach_1.5.2           generics_0.1.3          RCurl_1.98-1.12        
## [64] munsell_0.5.0           scales_1.2.1            globals_0.16.2         
## [67] class_7.3-22            glue_1.6.2              tools_4.4.0            
## [70] data.table_1.14.8       ModelMetrics_1.2.2.2    gower_1.0.1            
## [73] grid_4.4.0              missForest_1.5          tidyverse_2.0.0        
## [76] ipred_0.9-14            colorspace_2.1-0        nlme_3.1-163           
## [79] GenomeInfoDbData_1.2.11 cli_3.6.1               fansi_1.0.5            
## [82] S4Arrays_1.3.0          lava_1.7.2.1            dplyr_1.1.3            
## [85] pcaMethods_1.95.0       gtable_0.3.4            sass_0.4.7             
## [88] digest_0.6.33           SparseArray_1.3.0       htmltools_0.5.6.1      
## [91] lifecycle_1.0.3         hardhat_1.3.0           MASS_7.3-60.1

6 References

Appendix

Dekermanjian, J.P., Shaddox, E., Nandy, D. et al. Mechanism-aware imputation: a two-step approach in handling missing values in metabolomics. BMC Bioinformatics 23, 179 (2022). https://doi.org/10.1186/s12859-022-04659-1