bambinos / bambi

BAyesian Model-Building Interface (Bambi) in Python.
https://bambinos.github.io/bambi/
MIT License
1.08k stars 124 forks source link

interpret.plot_slopes ignores prob argument and az.rcParams['stats.hdi_prob' : ...] #748

Closed jt-lab closed 1 year ago

jt-lab commented 1 year ago

@GStechschulte, we noticed that plot_slope seems to ignore the hdi prob arguments. We noticed it due to mismatches between the plot and the tabulated results. Neither changing the prob argument not the az rcParams seems to have an effect:

fig, ax = bmb.interpret.plot_slopes(
    model=model,
    idata=idata,
    wrt="factor2",
    average_by=["factor3", "factor1"],
    subplot_kwargs={"main": "factor3", "panel": "factor1"},
    fig_kwargs={"figsize": (16, 4), "sharey": True},
    legend=False,
    eps = 1,
    prob = 0.90
)

fig, ax = bmb.interpret.plot_slopes(
    model=model,
    idata=idata,
    wrt="factor2",
    average_by=["factor3", "factor1"],
    subplot_kwargs={"main": "factor3", "panel": "factor1"},
    fig_kwargs={"figsize": (16, 4), "sharey": True},
    legend=False,
    eps = 1,
    prob = 0.99
)

import arviz as az
az.rcParams["stats.hdi_prob"] = 0.999

fig, ax = bmb.interpret.plot_slopes(
    model=model,
    idata=idata,
    wrt="factor2",
    average_by=["factor3", "factor1"],
    subplot_kwargs={"main": "factor3", "panel": "factor1"},
    fig_kwargs={"figsize": (16, 4), "sharey": True},
    legend=False,
    eps = 1
)

These plots should be different: image

Code to reproduce the observation

import pandas as pd
import bambi as bmb
from matplotlib.pylab import plt

df = pd.read_csv('simulated-data.csv')

model = bmb.Model('p(correct, count) ~ factor1*factor2*factor3 + (factor3+factor2|individual)',
                  df, family='binomial')
idata = model.fit(tune=500, draws=1000, target_accept=0.7, init='adapt_diag')

(and the the above plotting code)

Dataset: simulated-data.csv

Also found a little type in the doc string: az.rcParam["stats.hdi_prob"] --> az.rcParams["stats.hdi_prob"]

As always, many thanks in advance!

GStechschulte commented 1 year ago

@jt-lab Again, thanks for raising the issue and providing the code. Indeed, this is a bug. This is because the probs variable in def slopes() was never passed to the method .get_estimate() and thus probs defaulted to 0.94. After passing probs, I get the following plots:

fig, ax = bmb.interpret.plot_slopes(
    model=model,
    idata=idata,
    wrt="factor2",
    average_by=["factor3", "factor1"],
    subplot_kwargs={"main": "factor3", "panel": "factor1"},
    fig_kwargs={"figsize": (16, 4), "sharey": True},
    legend=False,
    eps = 1,
    prob = 0.90
)

fig, ax = bmb.interpret.plot_slopes(
    model=model,
    idata=idata,
    wrt="factor2",
    average_by=["factor3", "factor1"],
    subplot_kwargs={"main": "factor3", "panel": "factor1"},
    fig_kwargs={"figsize": (16, 4), "sharey": True},
    legend=False,
    eps = 1,
    prob = 0.99
)

import arviz as az
az.rcParams["stats.hdi_prob"] = 0.999

fig, ax = bmb.interpret.plot_slopes(
    model=model,
    idata=idata,
    wrt="factor2",
    average_by=["factor3", "factor1"],
    subplot_kwargs={"main": "factor3", "panel": "factor1"},
    fig_kwargs={"figsize": (16, 4), "sharey": True},
    legend=False,
    eps = 1
)

image image image

I will open a PR. Many thanks! :)

jt-lab commented 1 year ago

Many thanks for update & looking into this!

GStechschulte commented 1 year ago

Many thanks for update & looking into this!

PR has been opened!