Comparing conditional_effects from multiple models in one plot?

lakemonster · June 22, 2020, 2:07pm

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

Solomon · June 23, 2020, 10:06pm

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.

sofiaws · May 11, 2022, 8:50am

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.

Topic		Replies	Views
Question re brms: conditional_effects General brms	4	1016	July 4, 2022
Plotting model effects with brms Modeling mixed-model	8	5602	May 27, 2020
Brms::conditional_effects plots all 2-way interactions General conditional-effects	2	967	September 3, 2021
Difference between conditional_effects plot and predicted probability plot brms brms	3	3440	May 26, 2022
Conditional_effects: level over other fixed effects brms	1	872	March 9, 2020

Comparing conditional_effects from multiple models in one plot?

Related topics