Fitting multimodal data

Questions
1

I’m trying to fit the following distribution to which I only have access to the binned data:

image
image880×190 20.3 KB

I’ve tried the following mixture function:

def mixture_density(w, alpha, beta, scalling, x):
    logp = tt.log(w) + pm.Normal.dist(mu=alpha,sigma=beta).logp(x)
    return scalling * tt.sum(tt.exp(logp))

and the model I’ve been using is the following:

with pm.Model() as frequency_model:
    
    alpha_1 = pm.HalfNormal('Alpha_1', sigma= 1.)
    beta_1 = pm.HalfNormal('Beta_1', sigma= 1., )
    
    alpha_2 = pm.HalfNormal('Alpha_2', sigma= 1.)
    beta_2 = pm.HalfNormal('Beta_2', sigma= 1., )
    
    mix_1 = mixture_density(alpha_1, beta_1, x = frequency)
    mix_2 = mixture_density(alpha_2, beta_2, x = frequency)
    
    noise_1 = pm.HalfNormal('noise_1', sigma = 1)
    noise_2 = pm.HalfNormal('noise_2', sigma = 1)
        
    N_1 = pm.Normal('N_1', mu = mix_1, sigma=noise_1)
    N_2 = pm.Normal('N_2', mu = mix_2, sigma=noise_2)
    
    w = pm.Dirichlet('w', a=np.array([1.,1.]))
    
    likelihood = pm.Mixture('Likelihood',w=w, comp_dists=[N_1, N_2], observed=frequency_density)
           
    trace = pm.sample(draws=3000, chains = 4, tune=2000, target_accept=0.8)

The results from the nuts sampler have been inconclusive and the model can’t distinguish between both distributions. Is there a better way to do this? Or should I split the dataset and adjust each mode independently (I’ve had some success with that)
Thanks for your answers!!

1 Like