R/posthoc_Levelwise.R
posthoc_Levelwise.RdThis function is a wrapper based on emmeans, and needs a ordinary linear model produced by simple_model or a mixed effects model produced by mixed_model or mixed_model_slopes (or generated directly with lm, lme4 or lmerTest calls). It also needs to know the fixed factor(s), which should match those in the model and data table.
posthoc_Levelwise(Model, Fixed_Factor, P_Adj = "fdr", Factor, ...)a model object fit using simple_model or mixed_model or related.
one or more categorical variables, provided as a vector (see Examples), whose levels you wish to compare pairwise. Names of Fixed_Factor should match Fixed_Factor used to fit the model. When more than one factor is provided e.g. Fixed_factor = c("A", "B"), this function passes this on as specs = A|B (note the vertical | between the two Fixed_Factor) to emmeans to produce comparisons between each level A with each other listed separately at each level of B.
method for correcting P values for multiple comparisons. Default is set to false discovery rate ("fdr"), can be changed to "none", "tukey", "bonferroni", "sidak". See Interaction analysis in emmeans in the manual for emmeans.
old argument name for Fixed_Factor; retained for backward compatibility.
additional arguments for emmeans such as lmer.df or others. See help for sophisticated models in emmeans.
returns an "emm_list" object containing contrasts and emmeans through emmeans.
The function will generate level-wise comparisons (as described in Comparisons and contrasts in emmeans), i.e. comparison between of every level of one factor separately at each level of the other factor.
By default, P values are corrected by the FDR method (which can be changed). If the model was fit by transforming the quantitative response variable using "log", "logit", "sqrt" etc., results will still be on the original scale, i.e. type = "response" is the default; data will be back-transformed (check results to confirm this), and for log or logit see Transformations and link functions in emmeans, ratios will be compared.
The first part of the emmeans results has the estimated marginal means, SE and CI ($emmeans), which are generated from the fitted model, and not the original data table. The second part has the results of the comparisons ($contrasts).
#make a linear model first
CholMod <- mixed_model(data = data_cholesterol,
Y_value = "Cholesterol",
Fixed_Factor = c("Hospital", "Treatment"),
Random_Factor = "Subject")
#note quotes used only for fixed Fixed_Factor
#to get comparisons between different hospitals separately for each level of Treatment
posthoc_Levelwise(Model = CholMod,
Fixed_Factor = c("Hospital", "Treatment"))
#> $emmeans
#> Treatment = After_drug:
#> Hospital emmean SE df lower.CL upper.CL
#> Hosp_a 147 23.3 20.3 98.4 196
#> Hosp_b 142 23.3 20.3 93.8 191
#> Hosp_c 105 23.3 20.3 56.4 154
#> Hosp_d 163 23.3 20.3 114.2 211
#> Hosp_e 162 23.3 20.3 113.8 211
#>
#> Treatment = Before_drug:
#> Hospital emmean SE df lower.CL upper.CL
#> Hosp_a 173 23.3 20.3 124.6 222
#> Hosp_b 161 23.3 20.3 111.9 209
#> Hosp_c 124 23.3 20.3 75.6 173
#> Hosp_d 162 23.3 20.3 113.1 210
#> Hosp_e 163 23.3 20.3 114.2 211
#>
#> Degrees-of-freedom method: kenward-roger
#> Confidence level used: 0.95
#>
#> $contrasts
#> Treatment = After_drug:
#> contrast estimate SE df t.ratio p.value
#> Hosp_a - Hosp_b 4.587 33 20.3 0.139 0.9898
#> Hosp_a - Hosp_c 42.004 33 20.3 1.273 0.6749
#> Hosp_a - Hosp_d -15.770 33 20.3 -0.478 0.8062
#> Hosp_a - Hosp_e -15.432 33 20.3 -0.468 0.8062
#> Hosp_b - Hosp_c 37.417 33 20.3 1.134 0.6749
#> Hosp_b - Hosp_d -20.356 33 20.3 -0.617 0.8062
#> Hosp_b - Hosp_e -20.019 33 20.3 -0.607 0.8062
#> Hosp_c - Hosp_d -57.774 33 20.3 -1.751 0.4841
#> Hosp_c - Hosp_e -57.436 33 20.3 -1.741 0.4841
#> Hosp_d - Hosp_e 0.338 33 20.3 0.010 0.9919
#>
#> Treatment = Before_drug:
#> contrast estimate SE df t.ratio p.value
#> Hosp_a - Hosp_b 12.662 33 20.3 0.384 0.9746
#> Hosp_a - Hosp_c 49.005 33 20.3 1.485 0.7088
#> Hosp_a - Hosp_d 11.432 33 20.3 0.347 0.9746
#> Hosp_a - Hosp_e 10.370 33 20.3 0.314 0.9746
#> Hosp_b - Hosp_c 36.343 33 20.3 1.102 0.7088
#> Hosp_b - Hosp_d -1.230 33 20.3 -0.037 0.9746
#> Hosp_b - Hosp_e -2.291 33 20.3 -0.069 0.9746
#> Hosp_c - Hosp_d -37.573 33 20.3 -1.139 0.7088
#> Hosp_c - Hosp_e -38.635 33 20.3 -1.171 0.7088
#> Hosp_d - Hosp_e -1.062 33 20.3 -0.032 0.9746
#>
#> Degrees-of-freedom method: kenward-roger
#> P value adjustment: fdr method for 10 tests
#>
#get comparisons between treatments separately at each hospital
posthoc_Levelwise(Model = CholMod,
Fixed_Factor = c("Treatment", "Hospital"))
#> $emmeans
#> Hospital = Hosp_a:
#> Treatment emmean SE df lower.CL upper.CL
#> After_drug 147 23.3 20.3 98.4 196
#> Before_drug 173 23.3 20.3 124.6 222
#>
#> Hospital = Hosp_b:
#> Treatment emmean SE df lower.CL upper.CL
#> After_drug 142 23.3 20.3 93.8 191
#> Before_drug 161 23.3 20.3 111.9 209
#>
#> Hospital = Hosp_c:
#> Treatment emmean SE df lower.CL upper.CL
#> After_drug 105 23.3 20.3 56.4 154
#> Before_drug 124 23.3 20.3 75.6 173
#>
#> Hospital = Hosp_d:
#> Treatment emmean SE df lower.CL upper.CL
#> After_drug 163 23.3 20.3 114.2 211
#> Before_drug 162 23.3 20.3 113.1 210
#>
#> Hospital = Hosp_e:
#> Treatment emmean SE df lower.CL upper.CL
#> After_drug 162 23.3 20.3 113.8 211
#> Before_drug 163 23.3 20.3 114.2 211
#>
#> Degrees-of-freedom method: kenward-roger
#> Confidence level used: 0.95
#>
#> $contrasts
#> Hospital = Hosp_a:
#> contrast estimate SE df t.ratio p.value
#> After_drug - Before_drug -26.165 4.13 20 -6.341 <.0001
#>
#> Hospital = Hosp_b:
#> contrast estimate SE df t.ratio p.value
#> After_drug - Before_drug -18.090 4.13 20 -4.384 0.0003
#>
#> Hospital = Hosp_c:
#> contrast estimate SE df t.ratio p.value
#> After_drug - Before_drug -19.164 4.13 20 -4.645 0.0002
#>
#> Hospital = Hosp_d:
#> contrast estimate SE df t.ratio p.value
#> After_drug - Before_drug 1.037 4.13 20 0.251 0.8041
#>
#> Hospital = Hosp_e:
#> contrast estimate SE df t.ratio p.value
#> After_drug - Before_drug -0.362 4.13 20 -0.088 0.9309
#>
#> Degrees-of-freedom method: kenward-roger
#>