Dear PYMC community,

I am building a simulation of how people rate images.

I am using a 5 parameter logistic regression model to estimate response.

The model works fine.

The problem starts when I try to compare groups to see how they differ in-terms of the logistic regression parameters.

I currently use 2 likelihood functions, but then I cannot compare models.

I would like to create a single likelihood function but different priors to the different groups, with some shared hyper priors.

I tried using aesara switch, with an index for group on top of the index for subject but it won’t run. I think the logic behind it wasn’t that good.

I create the data using this method

```
def add_sub(trials = 100, min_r = 0, max_r = 4, noise = .1,
slope=10, a1 = 1, a2 = 20, c_mu = 2, d1 = 10, d2 = 10, g_mu = 1,
axis = False):
x = np.arange(min_r, max_r, max_r/trials)
a = np.random.beta(a1, a2, 1)[0] # close to 0 but not 0.
b = np.random.normal(slope, 2, 1)[0] # this is where the magic happens
c = np.random.normal(c_mu, .5, 1)[0] # kind of average
d = np.random.beta(d1, d2, 1)[0] # close to 0.5
g = np.random.normal(g_mu, .2, 1)[0] # almost symetrical
# the 5PL formula
y = d + (a-d)/(1+(x/c)**b)**g + np.random.normal(0, noise, trials)
# replace negative values with 0 (you cannot rate lower than 0)
y = [0 if i < 0 else i for i in y]
df = pd.DataFrame({'x': x, 'y': y})
if axis:
return(x)
return(df)
```

and then try to extract the data using this method:

```
def pl5_multi_groups(df1, df2, idx1, idx2, n_subs1, n_subs2 = 0):
if n_subs2 == 0:
n_subs2 = n_subs1
with pm.Model() as five_PL_groups:
# hyper prior per group
mu_b = pm.Normal('mu_b', 6, 3, shape = 2)
mu_c = pm.Normal('mu_c', 3.5, 1, shape = 2)
mu_g = pm.TruncatedNormal('mu_g', 1.1, .2, lower = .4, shape = 2)
a1 = pm.TruncatedNormal('a1', 1, .3, lower = 0, shape = 2)
a2 = pm.Normal('a2', 20, 3, shape = 2)
d1 = pm.TruncatedNormal('d1', 12, 2, lower = 0, shape = 2)
d2 = pm.TruncatedNormal('d2', 8, 2, lower = 0, shape = 2)
# priors
a_1 = pm.Beta('a_1', a1[0], a2[0], shape = n_subs1)
b_1 = pm.Normal('b_1', mu_b[0], 2, shape = n_subs1)
c_1 = pm.Normal('c_1', mu_c[0], 2, shape = n_subs1)
d_1 = pm.Beta('d_1', d1[0], d2[0], shape = n_subs1)
g_1 = pm.Normal('g_1', mu_g[0],.2, shape = n_subs1)
a_2 = pm.Beta('a_2', a1[1], a2[1], shape = n_subs2)
b_2 = pm.Normal('b_2', mu_b[1], 2, shape = n_subs2)
c_2 = pm.Normal('c_2', mu_c[1], 2, shape = n_subs2)
d_2 = pm.Beta('d_2', d1[1], d2[1], shape = n_subs2)
g_2 = pm.Normal('g_2', mu_g[1],.2, shape = n_subs2)
eps = pm.Exponential('eps', 5, shape = 2)
# model prediction using the 5PL regression
y_hat1 = d_1[idx1] + (a_1[idx1]-d_1[idx1])/(1+(df1['x'].values/c_1[idx1])**b_1[idx1])**g_1[idx1]
y_hat2 = d_2[idx2] + (a_2[idx2]-d_2[idx2])/(1+(df2['x'].values/c_2[idx2])**b_2[idx2])**g_2[idx2]
# Likelihood (sampling distribution) of observations
rating1 = pm.Normal('rating1',y_hat1, eps[0], observed=df1.y)
rating2 = pm.Normal('rating2',y_hat2, eps[1], observed=df2.y)
trace = jx.sample_numpyro_nuts(tune = 3000, draws = 1000, target_accept=0.95)
return trace
```

my question is how to send the function which group each participant is and then to have the function select the priors based on the group.

PYMC 4.1.7

Arviz 0.12.1

JAX 0.3.17

Aesara 2.8.2

complete code