Marginal distribution estimation

Thank you!
I try to something like (Obviously, that is just pseudocode):

with pm.Model() as model1:
teta1 = pm.Beta(‘teta1’,1, 2)
teta2 = pm.Beta(‘teta2’,1, 2)
teta3 = pm.Beta(‘teta3’,1, 2)
teta4 = pm.Beta(‘teta4’,1, 2)

teta100 = pm.Beta(‘teta100’,1, 2)

mnorm = pm.MvNormal('mnorm', mu=[theta1,theta2,...,theta100], cov=cov)


mnorm.sample(10000)

Where parameters for each Beta are obtained by experiment.
I’m not sure that is exists easy solution to get that analytically.