Issue with observed data in joint distribution

Thanks a lot for the quick response. I have changed what you have told me, but I am not sure how to pass the data using DensityDist. Since it requires a function, I have done the following:

data = [shift,end_hour,lengths,dows,hrs]
with pm.Model() as total_model:
sd1 = pm.Uniform(‘sd_SD’, 0.25, 5, shape=num_shifts)
lambdas1 = pm.Uniform(‘lambdas_SD’, 0, 1.0, shape=num_shifts)
p_sd = lambda d: pm.NormalMixture(‘P_SD’, mu=np.arange(num_shifts), w=lambdas1, sd=sd1, observed=d)

        sd2 = pm.Uniform('sd_ED', 0.25, 5, shape=num_shifts)
        lambdas2 = pm.Uniform('lambdas_ED', 0, 1.0, shape=num_shifts)
        p_ed = lambda d: pm.NormalMixture('P_ED', mu=np.arange(num_shifts), w=lambdas2, sd=sd2, observed=d)

        sd3 = pm.Uniform('sd_L', 0.25, 5, shape=num_lenghts)
        lambdas3 = pm.Uniform('lambdas_L', 0, 1.0, shape=num_lenghts)
        p_length = lambda d: pm.NormalMixture('P_length', mu=np.arange(num_lenghts), w=lambdas3, sd=sd3, observed=d)

        prior_dow = pm.Dirichlet('prior_dow', a=pm.floatX((1.0 / num_dow) * np.ones(num_dow)))
        p_dow = lambda d: pm.Multinomial('P_dow',num_dow, prior_dow, shape=num_dow, observed=d)

        mean_hr = pm.Uniform("mh", 50, 110)
        std_hr = pm.Uniform("sh", 1, 30)
        p_hr = lambda d: pm.Normal('P_hr',mean_hr, std_hr, observed=d) 

        pm.DensityDist('mylike', lambda d: p_sd(d[0]) * p_ed(d[1]) * p_length(d[2]) * p_dow(d[3]) * p_hr(d[4]), observed=data)
        trace = pm.sample(1000, tune=1000)

The error is:

ValueError: Input dimension mis-match. (input[0].shape[0] = 5, input[1].shape[0] = 7)

Respect to multiplying all the variables, I know the way they are described at this moment they are independent and having the distributions separately would be more efficient. In fact, I have them in this way before, but I want to add other factors that will change the independencies and thus I wanted to make it work in this way before doing other modifications.

Thanks