Gabriel Weindel writes:
I am a PhD student in psychology and I have a question about Bayesian statistics. I want to compare two posterior distributions of parameters estimated from a (hierarchical) cognitive model fitted on two dependent variables (hence both fits are completely separated). One fit is from a DV allegedly containing psychological process X and Y, and the other one is from a DV that only contains X. The test is to look up whether the cognitive model does ‘notice’ the removal of Y selectively in the parameter that is supposed to contain this process !
My take was to assume that, as I have access to the posterior distribution of the population parameters for both fits, I can simply compute the overlap (or equivalent) between both posterior distributions and if this overlap is high/low-to-null, conclude that there is high/low-to-no evidence that the true parameters of the fit on the two DVs are the same.
But my senior co-authors disagree with me and reviewers will probably also as, first this might be wrong and second this obviously goes against most of the statistics used in psychology and elsewhere were you need a criterion to decide between a null and an alternative hypothesis and where you rarely have access to a posterior distribution of the population parameter. However to me it appears to be both the most desirable and valid solution.
Does this reasoning seem valid to you ?
My quick answer is that I don’t think it makes sense to compare posterior distributions. Instead I think you should fit one larger model that includes both predictors.
Weindel responds:
I don’t see why it doesn’t make sense. We had thought about fitting a larger model but we would then add a dummy variable (DV1 = 0, DV2 = 1) and the two predictors would be highly correlated as they share a process (r = .85), wouldn’t that be a problem also?
My reply: Sure, when two predictors are highly correlated, then it’s hard from the data alone to tell them apart. That’s just the way it is!
> My quick answer is that I don’t think it makes sense to compare posterior distributions. Instead I think you should fit one larger model that includes both predictors.
Would Bayesian Workflow perspective suggest to fit all three models (A, B and A+B) and compare them all together?
I find it a little hard to wrap my head around the proposed inference without understanding exactly how the cognitive model represents psychological processes X and Y, but I will share my two cents anyway.
If there are two outcome variables which might or might not depend on the same underlying parameters, one simple way to test how well the shared paremeters between processes X and Y model the outcomes is to fit a restricted model (same set of parameters for X and Y) and compare it to a unrestricted model (different set of parameters for outcome X and Y) using, e.g., PSIS-LOO.
Let’s say that process X is parameterized by theta_x; and process Y is parameterized by theta_y. In the first, restricted, model, you might fit something like
DV_xy ~ f(theta_x, theta_y)
DV_x ~ g(theta_x)
In the unrestricted model, you might model each DV with a separate set of parameters, let’s say theta_x and theta_y for DV_xy and theta_x* for DV_x:
DV_xy ~ f(theta_x, theta_y)
DV_x ~ g(theta_x*)
If the restricted model is not considerably worse than the unrestricted model using LOO statistics, you have some evidence that the shared theta_x parameters is able to capture the underlying process to some degree.
I really don’t understand what’s being proposed, but I think it’s something along the lines of
y_a = f(a1) in dataset 1
vs
y_ab = f(a2,b) in dataset 2
and then compare the marginal distribution of a1 to a2
but what I absolutely think you should do is do:
y = f(a,b) (for the full dataset)
where b = 0 for those cases where it was supposedly set to zero by the experimenter.
anything else seems to make no sense.
I endorse Andrew’s suggestion. Although I’m not sure of the specific models the correspondent is interested in, I’ll try paraphrasing their scenario to see if I have it right:
There are two candidate models for some task, A and AB. Model AB contains all the same mechanisms (along with their corresponding parameters) as model A, plus something extra (“B”). It sounds like the question has to do with a parameter shared between the models, we’ll call it parameter V. Does the estimated value of V differ when B is included in the model (model AB) vs. when it is not (model A)?
Andrew’s suggestion, which I endorse, is to jointly fit both model A and model AB. For parameter V, you would then have a hierarchy that would look something like this:
1) Estimate of V for participant K conditional on model J
2) Distribution of V conditional on model J
3) Distribution of V
The research question can be addressed by directly subtracting the estimates at layer 2, i.e., (V conditional on model A) – (V conditional on model AB). The correlations in parameter estimates for V between the two models is actually the reason you want to do it this way. If you just compare the marginal distributions from separate model fits, they may look overlapping. But by looking directly at the difference distribution, you eliminate all the shared variance and focus specifically on the differences that depend on whether B is included in the model or not. For example, it might be the case that for every participant, the estimate of V under model A is higher than under model AB. But if there is considerable variability in posterior on V, the marginal difference may look small and will miss the consistent (but small) within-participant difference.