mami: Model averaging (and model selection) after multiple...

Description Usage Arguments Details Value Author(s) References See Also Examples

Description

Performs model selection/averaging on multiply imputed data and combines the resulting estimates. The package also provides access to less frequently used model averaging techniques and offers integrated bootstrap estimation.

Usage

1
2
3
4
5
6
7
8
9
mami(X, missing.data=c("imputed","none","CC"),  
  model=c("gaussian","binomial","poisson","cox"), outcome=1, id=NULL, 
  method=c("MA.criterion","MS.criterion","LASSO","LAE","MMA","JMA"), 
  criterion=c("AIC","BIC","BIC+","CV","GCV"), kfold=5, cvr=F,
  inference=c("standard","+boot"), B=20, X.org=NULL, CI=0.95,
  var.remove=NULL, aweights=NULL, add.factor=NULL, add.interaction=NULL, 
  add.transformation=NULL, add.stratum=NULL, ncores=1,    
  candidate.models=c("all","restricted","very restricted"), screen=0,
  report.exp=FALSE, print.time=FALSE, print.warnings=TRUE, ...)

Arguments

X

Either a list of multiply imputed datasets (each of them a dataframe), or an object of class ‘amelia’ created by Amelia II, or an object of class ‘mids’ created by mice, or a single dataframe.

missing.data

A character string, typically "imputed" when multiply imputed data are provided under X, or "CC" if a complete case analysis is desired, or "none" if there is no missing data.

model

A character string specifying the model family, either "gaussian" for linear regression models, "binomial" for logistic regression models, "poisson" for Poisson regression models, or "coxph" for Cox's proportional hazards models.

outcome

A character vector or integer specifying the variable(s) or column(s) which should be treated as outcome variable(s).

id

A character vector or integer specifying the variable or column to be used for a random intercept in the analysis model.

method

A character string specifying the model selection or model averaging technique: "MA.criterion" for model averaging based on exponential AIC or BIC weights, "MS.criterion" for stepwise variable selection based on AIC or BIC; or (all subset) model selection based on (generalized) cross validation, "MMA" for Mallow's model averaging and "JMA" for Jackknife Model Averaging (only linear model, see mma and jma), "LASSO" for model selection with LASSO and additional lasso averaging (see lae) if method="LAE" ("LAE" not available for Cox model).

criterion

A character string specifying the model selection criterion used for criterion-based model selection and averaging; currently either "AIC" for Akaike's Information Criterion, "BIC" for the Bayes criterion of Schwarz, "BIC+" for the Bayes criterion of Schwarz with quicker model averaging based on the leaps algorithm of the BMA package. Model selection can also be done based on the cross validation error (CV, based on the squared loss function), or GCV for generalized cross validation.

kfold

An integer specifying kfold cross validation; to be used when applying shrinkage estimation (method="LASSO", method="LAE") or criterion CV is used for model selection.

cvr

A logical value specifying whether the subsets used for cross validation should be shuffled or not.

inference

A character string, either "standard" for applying multiple imputation combining rules upon the model selection/averaging estimates or "+boot" if additional bootstrap based inference is required. See reference section for more information.

B

An integer indicating the number of bootstrap replications to use (when inference = "+boot" is chosen).

X.org

A dataframe consisting of the original unimputed data (which needs to be specified when inference = "+boot" is chosen).

CI

A scalar greater than 0 and less than 1 specifying the confidence of the confidence interval.

var.remove

Either a vector of character strings or integers, specifying the variables or columns which are part of the data but not to be considered in the model selection/averaging procedure.

aweights

A weight vector that is relevant to the analysis model.

add.factor

Either a vector of character strings or integers, specifying the variables which should be treated as categorical/factors in the analysis. Variables which are already defined to be factors in the data are detected automatically and do not necessarily need to be specified with this option.

add.interaction

A list of either character strings or integers, specifying the variables which should be added as interactions in the analysis model.

add.transformation

A vector of character strings, specifying transformations of variables which should be added to the analysis models.

add.stratum

A character vector or integer specifying the variable used as a stratum in Cox regression analysis.

ncores

An integer specifying the number of cores (threads) to be utilized.

candidate.models

A character string specifying whether for criterion based model selection/averaging all possible candidate models should be considered ("all"), or only candidate models with a limited amount of variables ("restricted","very restricted").

screen

Number of variables which should be removed in an initial screening step with LASSO.

report.exp

A logical value specifying whether exponentiated coefficients should be reported or not.

print.time

A logical value specifying whether analysis time and anticipated estimation time for bootstrap estimation should be printed.

print.warnings

A logical value specifying whether warnings and any other output from the function should be printed or not.

...

Further arguments to be passed, i.e. for functions lae or jma; dredge from the MuMIn package; bic.glm and bic.surv from the BMA package; or cv.glmnet from package glmnet.

Details

Model selection/averaging will be performed on each imputed dataset. The results will be combined according to formulae (7)-(10) in Schomaker and Heumann (CSDA, 2014), see References below for more details. If inference="+boot" is chosen, then the procedure described in Table 1 will be performed in addition to standard MI inference. For longitudinal data (specified via id) the bootstrap is based on the subject/person/id level. To obtain insightful results from bootstrap estimation B should be large, at least B>200 and plot.mami may be used.

Note that a variable will be formally selected if it is selected (by means of either model selection or averaging) in at least one imputed set of data, but its overall impact will depend on how often it is chosen. As a result, effects of variables which are not supported throughout imputed datasets and candidate models will simply be less pronounced. Variable importance measures based on model averaging weights are calculated for each imputed dataset and will then be averaged.

If method="MA.criterion" is chosen and the number of variables is large, then computation time might be a burden. The package provides several options to decrease computation time, see Section 6.1 of the manual at http://mami.r-forge.r-project.org/. One option is to parallize computation using the option ncores.

The function provides access to linear, logistic, Poisson and Cox proportional hazard models; one may add a random intercept to each of these models with the id option. Other models are not supported yet. Variables used for the imputation model but not needed for the analysis model can be removed with option var.remove.

Value

Returns an object of class ‘mami’:

coefficients.ma

A matrix of coefficients, standard errors and confidence intervals for model averaging estimators.

coefficients.ma.boot

A matrix of coefficients and bootstrap results (confidence intervals, mean, standard error) for model averaging.

coefficients.s

A matrix of coefficients, standard errors and confidence intervals for model selection estimators.

coefficients.boot.s

A matrix of coefficients and bootstrap results (confidence intervals, mean, standard error) for model selection.

variable.importance

A vector containing the variable importance for each variable based on model averaging weights.

boot.results

A list of detailed estimation results for each bootstrap sample. The first list element refers to the results from model selection, the second entry the results from model averaging.

Author(s)

Michael Schomaker

References

Schomaker, M., Heumann, C. (2014) Model Selection and Model Averaging after Multiple Imputation, Computational Statistics & Data Analysis, 71:758-770

See Also

plot.mami to visualize bootstrap results, lae and mma for model averaging techniques.

Examples

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
citation("MAMI")

####################################################
# Example 1: Freetrade example from Amelia package #
#            Cross-Section-Time-Series Data        #
#            Linear and linear mixed model         #
####################################################
set.seed(24121980)
library(Amelia)
data(freetrade)
freetrade$pop <- log(freetrade$pop) # in line with original publication
freetrade_imp <- amelia(freetrade, ts = "year", cs = "country", noms="signed", 
                        polytime = 2, intercs = TRUE, empri = 2)

# AIC based model averaging and model selection in a linear model after MI 
# (with and without bootstrapping)
mami(freetrade_imp, method="MS.criterion", outcome="tariff",add.factor=c("country"))
mami(freetrade_imp, method="MA.criterion", outcome="tariff",add.factor=c("country"))
m1 <- mami(freetrade_imp, method="MS.criterion", outcome="tariff",add.factor=c("country"),
           inference="+boot",B=25,X.org=freetrade,print.time=TRUE)

m1            # be patient with bootstrapping, increase B>=200 for better results
summary(m1)   # easier to read
plot(m1, plots.p.page="4")

# For comparison: Mallow's model averaging (MMA) and complete case analysis 
mami(freetrade_imp, method="MMA", outcome="tariff",add.factor=c("country"))
mami(freetrade, method="MS.criterion", missing.data="CC", outcome="tariff",
     add.factor=c("country")) #Note the difference to imputed analysis (e.g. usheg)     

# Use linear mixed model with random intercept for country using "id"
mami(freetrade_imp, outcome="tariff", id="country")

####################################################################
# Example 2: HIV treatment data, linear model and Cox model        #
####################################################################

# Impute with Amelia
data(HIV)
HIV_imp <- amelia(HIV, m=5, idvars="patient",noms=c("hospital","sex","dead","tb","cm"),
                  ords=c("period","stage"),logs=c("futime","cd4"),
                  bounds=matrix(c(3,7,9,11,0,0,0,0,3000,5000,200,150),ncol=3,nrow=4))

# i)  Cox PH model
# Model selection (with AIC) to select risk factors for the hazard of death, 
# reported as hazard ratios
# Also: add transformations and interaction terms to candidate models 
mami(HIV_imp, method="MS.criterion",model="cox", outcome=c("futime","dead"), 
     add.factor=c("hospital","stage","period"), add.transformation=c("cd4^2","age^2"), 
     add.interaction=list(c("cd4","age")), report.exp=TRUE, var.remove=c("patient","cd4slope6"))
      
# Similar as above (no hazard ratios reported, no interaction, hospitals as stratum),
# but with boostrap CI and visualization of bootstrap distribution 
m2 <- mami(HIV_imp, method="MS.criterion",model="cox", inference="+boot",
           X.org=HIV, outcome=c("futime","dead"), add.factor=c("stage","period"),
           add.transformation=c("cd4^2","age^2"), add.stratum="hospital", B=25, 
           var.remove=c("patient","cd4slope6"),report.exp=TRUE,
           print.time=TRUE,print.warnings=FALSE) 
summary(m2)
plot(m2)

# Model averaging
# 3 variables are initially removed via a LASSO screening step
m3 <- mami(HIV_imp, method="MA.criterion",model="cox", outcome=c("futime","dead"),
           add.factor=c("stage","period"),add.transformation=c("cd4^2","age^2"),
           add.stratum="hospital", var.remove=c("patient","cd4slope6"),screen=3) 
summary(m3)

# ii) Linear model
# Model selection and averaging to identify predictors for immune recovery 6 months
# after starting antiretroviral therapy, presented as CD4 slope which is the average
# change in number of CD4 cells per week (deaths are ignored for this example)

# AIC based model selection (stepAIC) after multiple imputation 
mami(HIV_imp, method="MS.criterion", outcome="cd4slope6", 
     add.factor=c("hospital","stage","period"), var.remove=c("patient","dead","futime"))


#########################################################################################
# Example 3:   Model selection/averaging with no missing data, using shrinkage          #
# Example from Tibshirani, R. (1996) Regression shrinkage and selection via the lasso,  #
# Journal of the Royal Statistical Society, Series B 58(1), 267-288.                    #
# Useful to use mami to obtain Bootstrap CI after model selection/averaging             #
#########################################################################################

data(Prostate)
mami(Prostate,method="LASSO",missing.data="none",outcome="lpsa" 
     ,kfold=10) # LASSO (selection) based on 10-fold CV
mami(Prostate,method="LAE",missing.data="none",outcome="lpsa" 
     ,kfold=10) # LASSO averaging based on 10-fold CV
m4 <- mami(Prostate,missing.data="none",outcome="lpsa", inference="+boot",
           B=50,print.time=TRUE) # AIC based averaging with Boostrap CI
summary(m4)
plot(m4) # a few bimodal distributions: effect or not?


###################################################
# Example 4: use utilities from other packages    #
###################################################

# Examples (without missing data for simplicity)
# Model Averaging with AIC: use restrictions as done in "dredge"
# Example: Candidate models cannot contain "svi" and "lcp" at the same time
mami(Prostate,missing.data="none",outcome="lpsa", subset = !(svi && lcp))   
# Example: Make Occam's Window smaller when using BMA package
mami(Prostate,missing.data="none",outcome="lpsa",criterion="BIC+", OR=5)  
# Example: Use Elastic Net instead of Lasso (passing on alpha to glmnet())
mami(Prostate,missing.data="none",outcome="lpsa",method="LASSO",alpha=0.5)

MAMI documentation built on May 6, 2019, 3:02 p.m.