Comparing conditional_effects from multiple models in one plot?

I am fitting dozens of models to the same set of population and group effects. I’m interested in determining the effect of genetic admixture on a set of traits. Instead of presenting the reader with dozens of conditional_effects() plots, I would like to present one conditional effect plot with all of the effects of ‘admixed’ in the same space to easily and tidily show the conditional effects of admixture.

I’ve thought of two solutions to this problem, but neither seems ‘right’ and I’d like to ask for some guidance from those on this forum.

The first solution is to fit models with rescaled responses, and to then plot the conditional effects in the same figure. Does this make sense? I’m rescaling responses by 2*sd as suggested by Gilman (2008).

Here’s an example to run

z.dat <- data.frame(
ind_code = c("USDA16_02","JKH_10","DPR_12","OFR_15","DCK_10","CPL_02","TUR_08","CHL_12","FNO_11","MSG_02","OFR_15","NEG_08","RMP_13","TUR_08","HST_01","HST_12","HBY_07","LPD_11","OFR_03","CYH_15","USDA12_07","OUT_03","OFR_06","HBY_07","CLK_01","TUR_05","MSG_01","LON_07","TIM_03","CPL_08","USDA10_07","UMI_03","DPR_01","CLK_03","KAP_03","FNO_15","USDA18_08","JKH_02","JKH_10","FNO_08","SSR_08","NEG_15","KAP_14","FNO_04","USDA13_04","SSR_01","NEG_08","LPD_13","NEG_10","SSR_12","HWK_05","USDA14_09","LPD_13","NBY_05","OFR_03","LON_10","NBY_10","TIM_14","NEG_11","TUR_01","TIM_05","MMT_08","MMT_10","USDA7_10","OUT_14","FNO_14","HBY_12","CPL_09","TIM_12","USDA14_09","RMP_11","DPR_15","MSG_15","SLC_04","CHL_13","CPL_12","USDA2_10","LPD_07","GAM_01","SLC_05","LON_07","MSG_02","USDA17_10","LPD_09","MMT_08","USDA7_09","JKH_03","USDA12_02","USDA12_10","OFR_09","DPR_09","USDA9_03","RMP_03","USDA2_01","CYH_15","CYH_10","OUT_11","USDA16_02","FNO_11","HWK_01"),
rep_code = c("USDA16_02_R2","JKH_10_R3","DPR_12_R3","OFR_15_R3","DCK_10_R32","CPL_02_R1","TUR_08_R4","CHL_12_R2","FNO_11_R1","MSG_02_R2","OFR_15_R2","NEG_08_R3","RMP_13_R2","TUR_08_R5","HST_01_R1","HST_12_R1","HBY_07_R2","LPD_11_R3","OFR_03_R2","CYH_15_R3","USDA12_07_R2","OUT_03_R5","OFR_06_R2","HBY_07_R2","CLK_01_R2","TUR_05_R5","MSG_01_R3","LON_07_R1","TIM_03_R2","CPL_08_R5","USDA10_07_R3","UMI_03_R5","DPR_01_R1","CLK_03_R2","KAP_03_R3","FNO_15_R2","USDA18_08_R2","JKH_02_R2","JKH_10_R1","FNO_08_R5","SSR_08_R3","NEG_15_R3","KAP_14_R4","FNO_04_R6","USDA13_04_R1","SSR_01_R3","NEG_08_R2","LPD_13_R3","NEG_10_R3","SSR_12_R3","HWK_05_R6","USDA14_09_R2","LPD_13_R2","NBY_05_R2","OFR_03_R3","LON_10_R3","NBY_10_R3","TIM_14_R3","NEG_11_R3","TUR_01_R5","TIM_05_R4","MMT_08_R22","MMT_10_R1","USDA7_10_R1","OUT_14_R4","FNO_14_R5","HBY_12_R2","CPL_09_R1","TIM_12_R2","USDA14_09_R2","RMP_11_R2","DPR_15_R1","MSG_15_R2","SLC_04_R1","CHL_13_R3","CPL_12_R2","USDA2_10_R3","LPD_07_R1","GAM_01_R2","SLC_05_R2","LON_07_R3","MSG_02_R1","USDA17_10_R2","LPD_09_R3","MMT_08_R1","USDA7_09_R2","JKH_03_R3","USDA12_02_R1","USDA12_10_R2","OFR_09_R2","DPR_09_R3","USDA9_03_R1","RMP_03_R3","USDA2_01_R1","CYH_15_R3","CYH_10_R4","OUT_11_R4","USDA16_02_R2","FNO_11_R5","HWK_01_R3"),
admixed	 = c("unadmixed","admixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","admixed","admixed","unadmixed","admixed","admixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","admixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","admixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","admixed","unadmixed","admixed","admixed","admixed","unadmixed","admixed","unadmixed","unadmixed","unadmixed","admixed","admixed","unadmixed","admixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","admixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","admixed","unadmixed","unadmixed","unadmixed","unadmixed","admixed","unadmixed","admixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","admixed","unadmixed","unadmixed","unadmixed","unadmixed","admixed","unadmixed","unadmixed","unadmixed","unadmixed","unadmixed","admixed","unadmixed","unadmixed","admixed","admixed","unadmixed","admixed","unadmixed"),
col      = c(10,14, 3,20,12,19,12,14,12,22,17,15,16,12, 5, 6,22,14, 7, 5,23, 3,18,14, 3,17,18,24,24,14, 4, 2,13, 1,19, 8, 5,10,16,15,10,24, 2, 4,20, 4,11,12,13,13,11,21,10, 3,22, 1, 6,11,11,12, 8,19, 2,12,14, 7, 2, 8,15,16,16,11, 4,16,11,20,17,15,22,19,22,19,15,21, 4, 4,14, 5,25,13,24, 3, 8,12,24,18,20,18, 5,25),
row      = c(16,24,18,33,27, 8, 9,19,30,16,16,28, 9, 1,39,27,26,31,33, 8,17,26,36,28, 1,22,29,28,11,16, 2,27,21,26,20,22,36, 2,30,24,12, 9, 2,20,28,39,16,14,26,28,20,24,33,36,28,28, 7,28,18,19,35, 9,37,28,39,38,14,16,15,38,40,32, 8, 4, 5,40,34,23,29,20,24,37, 6,36,16, 6,29,26, 6, 9, 1,29,17,12, 1,15,29,39,26, 2),
z.X1	 = c(-0.197179571, 0.633328136,-0.112861270,-0.623179875,-0.320074872,-0.899832077,-0.115065670, 0.043651131,-0.409353073, 1.178917141,-0.298581972, 0.179772832, 0.695051337,-0.428641573,-0.259453872,-0.061608970, 0.083881431,-0.104594770,-0.078141970, 0.116396332,-0.494222474,-0.530043974,-0.516817574, 0.006176331, 0.054122031,-0.504693374, 1.328265242, 0.094352331, 0.335734134,-0.125536570,-0.298581972,-0.275435772, 0.188039332,-0.021378670,-0.335505672, 0.905571539, 0.026015931, 0.458629435, 0.347858334, 0.855421438, 0.689540337, 0.265193333,-0.627588675,-0.383451373, 0.372106734, 1.222454041,-0.116167870, 0.734730537,-0.597278175, 0.069552831,-0.184504271,-0.242920872, 0.383128734,-0.610504575, 0.071206131, 0.462487135, 0.062939631, 1.301261342,-0.273231372,-0.628139775, 0.131827132,-0.403842073,-0.611055675,-0.018623170,-0.296928672,-0.635855175,-0.301337472,-0.735604276,-0.165215771,-0.247880772,-0.521777474,-0.230796671, 0.377066634, 0.082228131, 0.334080834,-0.069875470,-0.139314071,-0.328341372,-0.433050373,-0.041218270,-0.105145870, 0.326365433,-0.263862672, 0.187488232,-0.188361971,-0.045075970, 1.545949744, 0.031526931,-0.359202973,-0.530595074, 2.176408150,-0.360856273, 0.604670936,-0.515715374,-0.333301272, 0.323058833,-0.152540471, 0.018300531,-0.682698675,-0.265515972),
z.X2	 = c(-0.10383653,-0.06669064,-0.40897302,-0.38973318,-0.37832933,-0.13311310, 0.20530131, 0.06005221, 0.50536942, 0.21770550, 0.18979607,-0.53138107, 0.37565890, 0.20146668,-0.15548733,-0.07642727,-0.44671912, 0.23291064,-0.75125540, 0.38752958,-0.77996511,-0.23584783,-0.40263755,-0.26142314,-0.24488422, 0.25091673, 0.12937564, 0.07689123,-0.12077559, 0.17739187,-0.54882030, 0.40486878,-0.68529977,-0.55268827,-0.61607637, 0.48652972,-0.48029713,-0.89563755, 0.07965883, 0.77416029,-0.79677079, 0.98336435,-0.03287921,-0.16228963, 0.19656502,-0.71167536, 0.86409069, 0.27142366, 1.05218761,-0.60777356, 1.09626918, 0.09846519, 0.40776976,-0.35485473,-0.28366399, 0.11767168, 0.17975934,-0.45132067, 0.52577632, 0.12157300, 0.61907452,-0.35462131,-0.54498567,-0.73368280, 0.85755515, 0.26738897, 1.05962346,-0.68223207,-0.12421009, 0.08602765,-0.14351661,-0.55242152,-0.75659054,-0.11807468,-0.17282652, 0.20550138, 0.49049772,-0.12501036, 0.18672836, 0.16608805,-0.39443477,-0.18836511, 0.78713134,-0.34401773, 0.12774176,-0.12854489,-0.35842260,-0.21484073,-0.53161448,-0.23761509, 0.19703185, 1.03234757, 0.23474460,-0.48663261, 0.56665680,-0.44718594, 1.67519823, 0.10076597, 0.98453141, 0.78149610),
z.X3	 = c(-0.15499560, 0.10418985,-0.35212257, 0.20640384,-0.15499560, 0.88904722,-0.11849061,-0.29006408, 1.26869915,-0.48719104, 0.05308286,-0.07833511, 0.20640384, 0.20640384,-0.11849061,-0.26086008,-0.50909404,-0.11849061, 0.10418985, 0.20640384,-0.19150059, 0.15529684, 0.05308286,-0.38132656,-0.40688006, 0.32687031,-0.22435509,-0.32291857,-0.03452912,-0.11849061,-0.26086008, 0.10418985,-0.68066751,-0.48719104,-0.26086008,-0.15499560,-0.48719104,-0.40688006, 0.38892880, 0.10418985,-0.32291857, 0.53129828, 0.38892880, 0.00927687,-0.46163755,-0.66241501, 0.61160926,-0.62225952, 2.25068348, 0.05308286, 0.20640384,-0.03452912,-0.40688006, 0.20640384,-0.43608405,-0.53099704, 0.20640384,-0.69892001, 1.01316419, 0.45828829, 0.05308286, 0.61160926, 0.10418985,-0.57845353, 0.26481183, 0.53129828, 0.32687031, 0.20640384,-0.19150059,-0.03452912, 1.41836962,-0.32291857,-0.03452912,-0.38132656,-0.55655053,-0.15499560, 1.13728117,-0.03452912, 0.20640384,-0.19150059,-0.35212257,-0.22435509, 0.20640384,-0.53099704,-0.07833511,-0.29006408,-0.53099704,-0.38132656,-0.19150059, 0.15529684,-0.82668748, 0.38892880,-0.03452912, 0.00927687, 0.69557075,-0.19150059, 1.58629259,-0.26086008, 1.58629259, 0.20640384)
)

# Fit the models
fit.zX1 <- brm(z.X1 ~ admixed + row + col + (1|ind_code),
               data = z.dat,
               family = gaussian(),
               chains = 4L, seed = 2398,
               control = list(max_treedepth = 15))

fit.zX2 <- brm(z.X2 ~ admixed + row + col + (1|ind_code),
               data = z.dat,
               family = gaussian(),
               chains = 4L, seed = 2398,
               control = list(max_treedepth = 15))

fit.zX3 <- brm(z.X3 ~ admixed + row + col + (1|ind_code),
               data = z.dat,
               family = gaussian(),
               chains = 4L, seed = 2398,
               control = list(max_treedepth = 15))


# Output the conditional effects.
cond_mat.zX1 <- conditional_effects(fit.zX1)$admixed %>% mutate(X = "zX1") 
cond_mat.zX2 <- conditional_effects(fit.zX2)$admixed %>% mutate(X = "zX2") 
cond_mat.zX3 <- conditional_effects(fit.zX3)$admixed %>% mutate(X = "zX3") 

# Change all the different names of the group estimates to X_mean
names(cond_mat.zX1)[2] <- "X_mean"
names(cond_mat.zX2)[2] <- "X_mean"
names(cond_mat.zX3)[2] <- "X_mean"

# Row bind the data.
# This object is a dataframe of the conditional effect estimates
# for admixed and unadmixed for each trait.
cond_df <- rbind(
 cond_mat.zX1, 
 cond_mat.zX2,
 cond_mat.zX3)

ggplot(cond_df, aes(x = effect1__, y = estimate__, group = X)) + geom_point() + geom_line() 

That code makes a figure like this one,

Rplot

So, question is…what am I doing wrong here? Does this kind of analysis make statistical sense, or am I violating some assumption I’m ignorant of? How can I compare the effect of admixed/unadmixed on multiple Y’s from different models jointly?

Thanks all for the comments,

K

System Info:

  • Operating System: Pop!_os
  • brms Version: 2.12.0

If I’m following, you mentioned you had at least two solutions, but it appears you only detailed the first. At any rate, your approach of putting the response variables on the same scale and then showing them all in one plot by combining the data from multiple conditional_effects plots makes sense, to me. Once you find a way to express the uncertainty in the effects with uncertainty intervals and mark off the different models with color or something else, you’ll be golden.

Hi @lakemonster, were you able to create this plot? I’m also interested in comparing conditional_effects from different models in the same plot!
Thank you!

Best,
S.