Analyse de puissance pour les modèles à effets mixtes à plusieurs niveaux / linéaires longitudinaux.
Le but de powerlmm est d'aider à concevoir des études de traitement longitudinal (groupes parallèles), avec ou sans regroupement de niveau supérieur (par exemple, regroupés longitudinalement par des thérapeutes, des groupes ou un médecin) et des données manquantes. Les principales caractéristiques du package sont:
powerlmm peut être installé à partir de Cran et Github.
# CRAN, version 0.4.0
install.packages( " powerlmm " )
# GitHub, dev version
devtools :: install_github( " rpsychologist/powerlmm " )Ceci est un exemple de mise en place d'un modèle longitudinal à trois niveaux avec des pentes aléatoires au niveau du sujet et du cluster, avec différents modèles de données manquants par bras de traitement. Des entrées standardisées relatives sont utilisées, mais les valeurs de paramètres brutes non standardisés peuvent également être utilisées.
library( powerlmm )
d <- per_treatment( control = dropout_weibull( 0.3 , 2 ),
treatment = dropout_weibull( 0.2 , 2 ))
p <- study_parameters( n1 = 11 ,
n2 = 10 ,
n3 = 5 ,
icc_pre_subject = 0.5 ,
icc_pre_cluster = 0 ,
icc_slope = 0.05 ,
var_ratio = 0.02 ,
dropout = d ,
effect_size = cohend( - 0.8 ,
standardizer = " pretest_SD " ))
p
# >
# > Study setup (three-level)
# >
# > n1 = 11
# > n2 = 10 x 5 (treatment)
# > 10 x 5 (control)
# > n3 = 5 (treatment)
# > 5 (control)
# > 10 (total)
# > total_n = 50 (treatment)
# > 50 (control)
# > 100 (total)
# > dropout = 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10 (time)
# > 0, 0, 1, 3, 6, 9, 12, 16, 20, 25, 30 (%, control)
# > 0, 0, 1, 2, 4, 5, 8, 10, 13, 17, 20 (%, treatment)
# > icc_pre_subjects = 0.5
# > icc_pre_clusters = 0
# > icc_slope = 0.05
# > var_ratio = 0.02
# > effect_size = -0.8 (Cohen's d [SD: pretest_SD]) plot( p )get_power( p , df = " satterthwaite " )
# >
# > Power Analyis for Longitudinal Linear Mixed-Effects Models (three-level)
# > with missing data and unbalanced designs
# >
# > n1 = 11
# > n2 = 10 x 5 (treatment)
# > 10 x 5 (control)
# > n3 = 5 (treatment)
# > 5 (control)
# > 10 (total)
# > total_n = 50 (control)
# > 50 (treatment)
# > 100 (total)
# > dropout = 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10 (time)
# > 0, 0, 1, 3, 6, 9, 12, 16, 20, 25, 30 (%, control)
# > 0, 0, 1, 2, 4, 5, 8, 10, 13, 17, 20 (%, treatment)
# > icc_pre_subjects = 0.5
# > icc_pre_clusters = 0
# > icc_slope = 0.05
# > var_ratio = 0.02
# > effect_size = -0.8 (Cohen's d [SD: pretest_SD])
# > df = 7.971459
# > alpha = 0.05
# > power = 68%Les tailles de cluster inégales sont également prises en charge, les tailles de cluster peuvent être aléatoires (échantillonnées), soit la distribution marginale peut être spécifiée.
p <- study_parameters( n1 = 11 ,
n2 = unequal_clusters( 2 , 3 , 5 , 20 ),
icc_pre_subject = 0.5 ,
icc_pre_cluster = 0 ,
icc_slope = 0.05 ,
var_ratio = 0.02 ,
effect_size = cohend( - 0.8 ,
standardizer = " pretest_SD " ))
get_power( p )
# >
# > Power Analyis for Longitudinal Linear Mixed-Effects Models (three-level)
# > with missing data and unbalanced designs
# >
# > n1 = 11
# > n2 = 2, 3, 5, 20 (treatment)
# > 2, 3, 5, 20 (control)
# > n3 = 4 (treatment)
# > 4 (control)
# > 8 (total)
# > total_n = 30 (control)
# > 30 (treatment)
# > 60 (total)
# > dropout = No missing data
# > icc_pre_subjects = 0.5
# > icc_pre_clusters = 0
# > icc_slope = 0.05
# > var_ratio = 0.02
# > effect_size = -0.8 (Cohen's d [SD: pretest_SD])
# > df = 6
# > alpha = 0.05
# > power = 44%Les tailles de cluster suivent une distribution de Poisson
p <- study_parameters( n1 = 11 ,
n2 = unequal_clusters( func = rpois( 5 , 5 )), # sample from Poisson
icc_pre_subject = 0.5 ,
icc_pre_cluster = 0 ,
icc_slope = 0.05 ,
var_ratio = 0.02 ,
effect_size = cohend( - 0.8 ,
standardizer = " pretest_SD " ))
get_power( p , R = 100 , progress = FALSE ) # expected power by averaging over R realizations
# >
# > Power Analyis for Longitudinal Linear Mixed-Effects Models (three-level)
# > with missing data and unbalanced designs
# >
# > n1 = 11
# > n2 = rpois(5, 5) (treatment)
# > - (control)
# > n3 = 5 (treatment)
# > 5 (control)
# > 10 (total)
# > total_n = 24.96 (control)
# > 24.96 (treatment)
# > 49.92 (total)
# > dropout = No missing data
# > icc_pre_subjects = 0.5
# > icc_pre_clusters = 0
# > icc_slope = 0.05
# > var_ratio = 0.02
# > effect_size = -0.8 (Cohen's d [SD: pretest_SD])
# > df = 8
# > alpha = 0.05
# > power = 48% (MCSE: 1%)
# >
# > NOTE: n2 is randomly sampled. Values are the mean from R = 100 realizations.Plusieurs fonctions de commodité sont également incluses, par exemple pour créer des courbes de puissance.
x <- get_power_table( p ,
n2 = 5 : 10 ,
n3 = c( 4 , 8 , 12 ),
effect_size = cohend(c( 0.5 , 0.8 ), standardizer = " pretest_SD " ))plot( x ) Le package comprend une méthode de simulation flexible qui facilite l'étude des performances de différents modèles. Par exemple, comparons la différence de puissance entre le LMM à 2 niveaux avec 11 mesures répétées, à faire un ancova au post-test. En utilisant sim_formula différents modèles peuvent être ajustés au même ensemble de données pendant la simulation.
p <- study_parameters( n1 = 11 ,
n2 = 40 ,
icc_pre_subject = 0.5 ,
cor_subject = - 0.4 ,
var_ratio = 0.02 ,
effect_size = cohend( - 0.8 ,
standardizer = " pretest_SD " ))
# 2-lvl LMM
f0 <- sim_formula( " y ~ time + time:treatment + (1 + time | subject) " )
# ANCOVA, formulas with no random effects are with using lm()
f1 <- sim_formula( " y ~ treatment + pretest " ,
data_transform = transform_to_posttest ,
test = " treatment " )
f <- sim_formula_compare( " LMM " = f0 ,
" ANCOVA " = f1 )
res <- simulate( p ,
nsim = 2000 ,
formula = f ,
cores = parallel :: detectCores( logical = FALSE )) Nous résumons ensuite les résultats en utilisant summary . Examinons spécifiquement les effets du traitement.
summary( res , para = list ( " LMM " = " time:treatment " ,
" ANCOVA " = " treatment " ))
# > Model: summary
# >
# > Fixed effects: 'time:treatment', 'treatment'
# >
# > model M_est theta M_se SD_est Power Power_bw Power_satt
# > LMM -1.1 -1.1 0.32 0.32 0.93 0.93 .
# > ANCOVA -11.0 0.0 3.70 3.70 0.85 0.84 0.84
# > ---
# > Number of simulations: 2000 | alpha: 0.05
# > Time points (n1): 11
# > Subjects per cluster (n2 x n3): 40 (treatment)
# > 40 (control)
# > Total number of subjects: 80
# > ---
# > At least one of the models applied a data transformation during simulation,
# > summaries that depend on the true parameter values will no longer be correct,
# > see 'help(summary.plcp_sim)'Nous pouvons également consulter un modèle spécifique, voici les résultats du LMM à 2 lvl.
summary( res , model = " LMM " )
# > Model: LMM
# >
# > Random effects
# >
# > parameter M_est theta est_rel_bias prop_zero is_NA
# > subject_intercept 99.00 100.0 -0.00560 0 0
# > subject_slope 2.00 2.0 -0.00310 0 0
# > error 100.00 100.0 0.00086 0 0
# > cor_subject -0.39 -0.4 -0.03200 0 0
# >
# > Fixed effects
# >
# > parameter M_est theta M_se SD_est Power Power_bw Power_satt
# > (Intercept) 0.0048 0.0 1.30 1.30 0.050 . .
# > time 0.0011 0.0 0.25 0.25 0.054 . .
# > time:treatment -1.1000 -1.1 0.32 0.32 0.930 0.93 .
# > ---
# > Number of simulations: 2000 | alpha: 0.05
# > Time points (n1): 11
# > Subjects per cluster (n2 x n3): 40 (treatment)
# > 40 (control)
# > Total number of subjects: 80
# > ---
# > At least one of the models applied a data transformation during simulation,
# > summaries that depend on the true parameter values will no longer be correct,
# > see 'help(summary.plcp_sim)' La fonctionnalité de base du package est également implémentée dans une application Web brillante, destinée aux utilisateurs qui connaissent moins R. Lancez l'application en tapant
library( powerlmm )
shiny_powerlmm()