Added rope to forestplot.py in forestplot()

[GWeindel]

The new argument rope adds a shaded region of practical equivalence for all displayed variables in forestplot, shade is controlled by rope_alpha (default =.5). Contrary to posteriorplot rope can only be a tuple of length 2.
Not big but it is highly useful for me hence I assume it should be useful for others too.

[ahartikainen]

Sounds good. Can you post example?

[GWeindel]

@ahartikainen Sure, by the way it allowed me to spot an error in the ymax definition.
I post a reproducible example and link a plot I use.

Reproducible example :

import pystan
import arviz as az

data = dict(N1=10, y1=2, N2=10, y2=6)

model = """//  Comparison of two groups with Binomial
data {
  int<lower=0> N1;
  int<lower=0> y1;
  int<lower=0> N2;
  int<lower=0> y2;
}
parameters {
  real<lower=0,upper=1> theta1;
  real<lower=0,upper=1> theta2;
}
model {
  theta1 ~ beta(1,1);
  theta2 ~ beta(1,1);
  y1 ~ binomial(N1,theta1);
  y2 ~ binomial(N2,theta2);
}
generated quantities {
  real oddsratio;
  oddsratio = (theta2/(1-theta2))/(theta1/(1-theta1));
}
"""

sm = pystan.StanModel(model_code=model)
fit = sm.sampling(data=data, iter=1000, chains=2)

az.plot_forest(fit, var_names=["theta1","theta2"], credible_interval=0.95, combined=True, rope=(-.1,.1),
               figsize=(7,7))

[aloctavodia]

I like the idea, but it seems to me that the user should be allowed to set a rope per variable.

[canyon289]

Without having reviewed the details of the code, I would vote for consistency between the various plots

…

On Tuesday, December 11, 2018, Osvaldo Martin ***@***.***> wrote: I like the idea, but it seems to me that the user should be allowed to set a rope per variable. — You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub <#448 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/AG4S4YOVTOp48XuKAEgD4D-QFIFtEkkSks5u398ngaJpZM4ZNCYm> .

[GWeindel]

Sure I get it, I have no use for it but some might and consistency is defintely important. I will adapt the rope code of posteriorplot() to forestplot() when I have a bit time and perhaps generalize it to the ridgeplot option. Thanks for this highly useful package by the way.

[aloctavodia]

Great! And thank you for this contribution.

[GWeindel]

I severely rewrote my initial commit. Now rope has the same behavior as posteriorplot. I provide a reproductible example and a plot I could use.

Reproductible example :

import pystan
import arviz as az
import stan_utility
import matplotlib.pyplot as plt

data = dict(N1=10, y1=2, N2=10, y2=6)

model = """//  Comparison of two groups with Binomial
data {
  int<lower=0> N1;
  int<lower=0> y1;
  int<lower=0> N2;
  int<lower=0> y2;
}
parameters {
  real<lower=0,upper=1> theta1;
  real<lower=0,upper=1> theta2;
}
model {
  theta1 ~ beta(1,1);
  theta2 ~ beta(1,1);
  y1 ~ binomial(N1,theta1);
  y2 ~ binomial(N2,theta2);
}
generated quantities {
  real oddsratio;
  oddsratio = (theta2/(1-theta2))/(theta1/(1-theta1));
}
"""

sm = pystan.StanModel(model_code=model)
fit = sm.sampling(data=data, iter=1000, chains=2)

rope = {'theta1' : [{'rope': (-.1, .1)}],
 'theta2' : [{'rope': (.2, .5)}],
 'oddsratio' : [{'rope': (10, 20)}]}

ropeList = (-.1, .1)

az.plot_forest(fit, var_names=["theta1","theta2"],credible_interval=0.95, combined=True, rope=rope, rope_values=True)
plt.show()
az.plot_forest(fit, var_names=["theta1","theta2"], credible_interval=0.95, combined=True, rope=ropeList,  figsize=(7,7))
plt.show()

[GWeindel]

However, I don't master classes so I might have introduced a lot of redundancy...

[GWeindel]

Faster reproducible and more complete example taken from the gallery :

size = (1, 50)

data = {'normal': np.random.randn(*size),
    'gumbel': np.random.gumbel(size=size),
    'student t': np.random.standard_t(df=6, size=size),
    'exponential': np.random.exponential(size=size)}

rope = {'normal': [{"rope":(-.5,.5)}],
        'gumbel': [{"rope":(0,.5)}],
        'student t': [{"rope":(-.5,1)}],
        'exponential': [{"rope":(-.5,1)}]
       }

az.plot_forest(data, rope=(-1,1))
plt.show()
az.plot_forest(data, rope=rope)
plt.show()
az.plot_forest(data, rope=rope, rope_values=False)
plt.show()

[aloctavodia]

I like it. I think the rope should not be labeled here. I know it is for the posterior plot, by for that plot the hpd is also labeled.

[GWeindel]

You're right, I was hesitating whether to keep it.

[canyon289]

When running `./scripts/lint.sh/ I'm getting the following error

• SRC_DIR=/home/tbayes/repos/arviz
• echo 'Checking documentation...'
  Checking documentation...
• python -m pydocstyle --convention=numpy /home/tbayes/repos/arviz/arviz/
  /home/tbayes/repos/arviz/arviz/plots/forestplot.py:246 in public method display_multiple_ropes:
  D102: Missing docstring in public method

[canyon289]

@GWeindel Thank you for the quick edits. I'll do another round of review either today or tomorrow.

[canyon289]

@GWeindel Thanks for the fixes. Unfortunately our CI checks are broken at this commit so I'm running things manually as you make changes. Here's the errors I'm getting.

If you'd like you can rebase your changes on top of the current arviz master to get CI checks to work again. If you need help with this let me know.

Alternatively you can run

./scripts/lint.sh

and it will generate the output below

(arviz) tbayes@Users-MBP-2:~/repos/arviz$ ./scripts/lint.sh
++ pwd

• SRC_DIR=/Users/tbayes/repos/arviz
• echo 'Checking documentation...'
  Checking documentation...
• python -m pydocstyle --convention=numpy /Users/tbayes/repos/arviz/arviz/
  /Users/tbayes/repos/arviz/arviz/plots/forestplot.py:246 in public method display_multiple_ropes:
  D210: No whitespaces allowed surrounding docstring text
  /Users/tbayes/repos/arviz/arviz/plots/forestplot.py:246 in public method display_multiple_ropes:
  D400: First line should end with a period (not 'd')
  /Users/tbayes/repos/arviz/arviz/plots/forestplot.py:246 in public method display_multiple_ropes:
  D401: First line should be in imperative mood ('Display', not 'Displays')

[ahartikainen]

LGTM.

I fixed some annoying pylint errors.

[canyon289]

Does not LGTM yet,
Need test for new multi rope functionality. I can add one later today

https://coveralls.io/builds/20848553/source?filename=arviz%2Fplots%2Fforestplot.py#L257
https://coveralls.io/builds/20848553/source?filename=arviz%2Fplots%2Fforestplot.py#L352

[canyon289]

Give me a couple more days, got caught up with other tasks

[canyon289]

Trying to add testing but messing up on arguments somewhere. The following code generates an exception

import arviz as az
import matplotlib.pyplot as plt

az.style.use('arviz-darkgrid')

centered_data = az.load_arviz_data('centered_eight')

rope = {'mu': [{'rope': (-.1, .1)}], 'tau': [{'rope': (.2, .5)}]}
_, axes = az.plot_forest(centered_data,
                         var_names=['mu', 'tau'],
                         rope=rope)
axes[0].set_title('Estimated theta for eight schools model')
plt.show()

will continue working on this

[canyon289]

This implementation seems to break when each variable has multiple chains

Combined = True works

import arviz as az
import matplotlib.pyplot as plt
az.style.use('arviz-darkgrid')
centered_data = az.load_arviz_data('centered_eight')
rope = {'mu': [{'rope': (-.1, .1)}], 'tau': [{'rope': (.2, .5)}]}
_, axes = az.plot_forest(centered_data,
                         var_names=['mu', 'tau'],
                         combined=True,
                         rope=rope)
axes[0].set_title('Estimated theta for eight schools model')
plt.show()

Combined = False fails

import arviz as az
import matplotlib.pyplot as plt
az.style.use('arviz-darkgrid')
centered_data = az.load_arviz_data('centered_eight')
rope = {'mu': [{'rope': (-.1, .1)}], 'tau': [{'rope': (.2, .5)}]}
_, axes = az.plot_forest(centered_data,
                         var_names=['mu', 'tau'],
                         combined=False,
                         rope=rope)
axes[0].set_title('Estimated theta for eight schools model')
plt.show()

I have a couple suggested paths but happy to hear any others

1. Figure out how to make solution work with multiple chains
2. Provide a warning that per var rope won't work if multiple chains

I won't have time to code solution 1 if we choose to go that route, but can do solution 2.
Depending on which version we go for I can write the appropriate tests

[GWeindel]

I won't have the time (and perhaps sufficient skills) for solution 1. Hence I would suggest providing a warning.

[ahartikainen]

So what fails is

... label[ticks == y][0]

Which would mean that we don't have label for all y. Maybe test if y is found?

edit.

Line 331

self.display_multiple_ropes(rope, ax, y, linewidth, label[ticks == y][0])

to

if (ticks == y).any():
    self.display_multiple_ropes(rope, ax, y, linewidth, label[ticks == y][0]) 

[canyon289]

What happens is that ArviZ is trying to match the label with the line. When the chains are compressed label position matches the y position in the numpy array so one value is found, we get it out of the array in index zero.

In multi chain the ticks, and the y position, are off every so slightly, .175 if I remember correctly. The equality operator doesn't match anything and the array is empty, hence the index error.

I hope this not so great description yields some insight. I can provide more detail if needed

[ahartikainen]

Ok

We need y_group variable which is the mean of y locations per parameter (or value to be plotted). VarHandler probably needs to calculate these before hand, because all the yielding is disabling to do that later?

[ahartikainen]

Combined True

Combined False

[canyon289]

Now we just need to write a couple more tests and we're good!

Additional changes made afterwards

[canyon289]

Thanks for the contribution @GWeindel

[GWeindel]

Great ! Thanks for the help with the multilple chains and the tests !