Stretched beta prior

Questions
1

I’m trying to implement a stretched beta prior so I can specify an informative prior over a correlation coefficient. This has been done elsewhere (such as in JASP, see pic) by using Beta(2,2) then applying a transform to change the bounds from (0,1) to (-1, 1). This is easy enough to implement manually, the only problem is that the posterior still represents the updated beta spanning (0,1). In my context it’s a bit of an annoyance to have to remember to do this transformation every time I access this variable. Is there an easy way to implement this, or might this be worth a request on GitHub for a StretchedBeta prior?
19

2

Just access the transformed variable instead. If you have

beta_pr = pm.Beta('beta_pr', 2, 2)
corr = pm.Determinstic('corr', 2 * beta_pr - 1)

Then corr will show up in the trace.

1 Like
3

Thanks, that was rather simple.

So this does work, although in my model (minimum working example below) it throws an error when sampling from the prior.

import pymc3 as pm
import numpy as np

def _fisher_transformation(r):
    '''Apply the Fisher transformation. Converts from r to z.'''
    return 0.5*np.log((1+r)/(1-r))

# make observed data
n_studies = 10
r_obs = np.random.uniform(low=-1, high=1, size=n_studies)
n_obs = np.random.randint(low=10, high=50, size=n_studies)
z_obs = _fisher_transformation(r_obs)

with pm.Model() as model:
    beta_pr = pm.Beta('beta_pr', 2, 2)
    r_population = pm.Deterministic('r_population', 2 * beta_pr - 1)
    
    study_sigma = pm.Uniform('study_sigma', 0, 10) 
    
    z_study = pm.Normal('z_study',
                            mu=_fisher_transformation(r_population),
                            sd=study_sigma, shape=n_studies)

    z_obs = pm.Normal('z_obs', mu=z_study, 
                          sd=1/pm.math.sqrt(n_obs-3),
                          observed=z_obs)
    
# works
with model:
    posterior = pm.sample()
    
# does not work
with model:
    posterior = pm.sample_prior_predictive()

with the error ValueError: Cannot resolve inputs for ['r_population']

It does work if you do a hack where you sample from the prior before you define the rest of the model. But it’s not clear why you can sample from the posterior (ie the model works), but you can’t sample from the prior.

# define prior only
with pm.Model() as model:
    r_pop = pm.Beta('r_pop', 2, 2)
    r_population = pm.Deterministic('r_population', (r_pop*2)-1)

# sample from prior
with model:
    prior = pm.sample_prior_predictive()

# define rest of model
with model:
    ...
4

Did you upgrade to 3.7? There are a few fix introduced in the new version that should make prior predictive much smoother.

1 Like
5

Yes! This now works. All the answers to my questions are simple. Ben must try harder / take a holiday :roll_eyes:

6

+1 to holiday :grin:

1 Like