How can we build a mixture of mixtures?

Questions
19

I agree that my mixture should be one dimensional (and that’s what I want).

But we have to admit that this experiment :

mu_g = np.arange(5,data_max,5) # mu_g size is 40 
lp = np.exp(ln_mix.logp(np.asarray([mu_g]))) # a 1 dimensional mixture build on 5 components

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

tells that the logp function is expecting a 5 dimension vector… no ?
And providing a 5 dimension vector suppress this error, but is it really computing the correct lopg value of the mixture? I don’t know, since this input format is also accepted :

array([[  5],
       [ 10],
       [ 15],
       [ 20],
       [ 25],
       [ 30],
       [ 35],
       ...

Which is not 5-dim, so I don’t know why it says that it expected a 5 dimension input vector !

And anyway I always have an error when trying to use the logp output as 1-D vector
(either in as the weights of a new mixture, or simply trying to access a given element of that output vector).

But, is the lopg function supposed to be used with a vector as an input is the first place ??
(expecting to produce a vector output).

I seems to be true for basic distributions :

ln = pm.Lognormal.dist(mu=1,sd=0.5)
ln.logp([12,13]).eval()
=> array([-7.12059352, -7.68887369])

But not with mixtures (!):

 y = pm.Mixture.dist(w=[0.5,0.5],comp_dists=pm.Normal.dist(mu=np.asarray([1,2]),sd=np.asarray([1,1]),shape=2))
y.logp(np.asarray([12])).eval()
=> array(-51.612058177694394)
y.logp(np.asarray([13])).eval()
=> array(-62.11207558372233)
 y.logp(np.asarray([12,13])).eval()
=>array(-113.72413376141672)

The logp function seems to apply a sum operator, which sounds weird since a logp_sum function also exists :
y.logp_sum(np.asarray([12,13])).eval()
=> array(-113.72413376141672)

This explains the lack of dimension of the lopg output!

But I still don’t know how to call the logp function on an array without an explicit Python loop that totally ruins the performance.