I’d been studying the concepts of probabilistic programming using Pymc3 in the past couple of week and I have a question on how to implement a particular model, I will describe it next:
I have a set of N data elements, and a matrix NxN defining the interactions between each other of those elements (the matrix is symmetric). I also have a property E for each element, I want to model this property as it follows:
What the model says is that, the energy of a element M is a weighted sum over all the other elements in the data set. The function d(M, M_i) is known ( is the NxN matrix I described before)
I want to model the weights alpha_i and sigma using Bayesian inference.
My try so far goes like this:
indx = np.arange(0, N ) # The index for the N elements def estimator( alpha ,x , sigma , ): # My idea is that this estimator returns the estimated value for the energy for a given set alpha, x , and sigma return T.sum( alpha*T.exp( -1/(2.0*sigma**2 )*x**2 ) ) with pm.Model() as hi_model: # Positive prior on sigma, N positive priors for the alphas sigma = pm.HalfCauchy( 'sigma' , beta = 1 , shape = 1 ) alpha1 = pm.HalfCauchy( 'alphas' , beta = 1 , shape = N ) # An error parameter sigma_e = pm.HalfCauchy('sigma_e' , beta = 1 ) tau_e = sigma_e**-2 # The energy estimator build using the function above. # D2s is an NxN matrix, whose column i is a vector of length N #representing the interaction of the i element with all the other N elements energy_est = estimator( alpha1 , D2s[indx] , sigma ) ee = pm.Deterministic( 'muu' , energy_est ) y_like = pm.Normal('y_like', mu= ee , tau=tau_e, observed= energy_smp[:N] )
Pymc3 allow me to build and sample from this model without errors, but I’m not being able to use the sampled alphas an sigma to “predict” the N energies I’m giving as observed, I’m off by many orders of magnitude. So I’m concerned if the model is being understood as I want or if I’m missing something and thats not the way to build this kind of model.
Thanks in advance for the help