<jupyter_start><jupyter_text>A look at model comparisonsee : https://docs.pymc.io/notebooks/Bayes_factor.htmlhttps://docs.pymc.io/notebooks/model_comparison.html<jupyter_code>from collections import OrderedDict
import numpy as np
import pymc3 as pm
import seaborn as sns
sns.set(font_scale=1.5)
import scipy.stats
import matplotlib.pyplot as plt

%matplotlib inline
data = scipy.stats.nbinom.rvs(50, 0.04, size=20)
with pm.Model() as modelp:
    # fit a poisson and negbin
    p_mean = pm.Uniform('p_mean', 1, 10000)
    p = pm.Poisson('p', mu=p_mean, observed=data)
    tracep = pm.sample(10000, tune=4000)

with pm.Model() as modelnb:
    # fit a poisson and negbin
    nb_mean = pm.Uniform('nb_mean', 1, 10000)
    nb_alpha = pm.Uniform('nb_alpha', 1, 1000)
    nb = pm.NegativeBinomial('nb', mu=nb_mean, alpha=nb_alpha, observed=data)
    tracenb = pm.sample(10000, tune=4000)
    
    
pm.traceplot(tracep, combined=True)
pm.traceplot(tracenb, combined=True)
ppcp = pm.sample_ppc(tracep, samples=1000, model=modelp)
ppcnb = pm.sample_ppc(tracenb, samples=1000, model=modelnb)

sns.distplot(ppcp['p'], label='p')
sns.distplot(ppcnb['nb'], label='nb')
sns.distplot(data, label='data')
plt.legend()
pm.stats.waic(tracep, model=modelp, progressbar=True)
pm.stats.waic(tracenb, model=modelnb, progressbar=True)
pm.stats.loo(tracep, model=modelp, progressbar=True)
pm.stats.loo(tracenb, model=modelnb, progressbar=True)
modelp.name = 'Poisson'
modelnb.name = 'NegBin'

mods = OrderedDict()
mods[modelp]= tracep
mods[modelnb]= tracenb

comp = pm.stats.compare(mods)
comp
pm.plots.compareplot(comp)<jupyter_output><empty_output><jupyter_text>And try other data<jupyter_code>data = scipy.stats.poisson.rvs(50, size=20)
with pm.Model() as modelp:
    # fit a poisson and negbin
    p_mean = pm.Uniform('p_mean', 1, 10000)
    p = pm.Poisson('p', mu=p_mean, observed=data)
    tracep = pm.sample(10000, tune=4000)

with pm.Model() as modelnb:
    # fit a poisson and negbin
    nb_mean = pm.Uniform('nb_mean', 1, 10000)
    nb_alpha = pm.Uniform('nb_alpha', 1, 1000)
    nb = pm.NegativeBinomial('nb', mu=nb_mean, alpha=nb_alpha, observed=data)
    tracenb = pm.sample(10000, tune=4000)
    
    
ppcp = pm.sample_ppc(tracep, samples=1000, model=modelp)
ppcnb = pm.sample_ppc(tracenb, samples=1000, model=modelnb)

sns.distplot(ppcp['p'], label='p')
sns.distplot(ppcnb['nb'], label='nb')
sns.distplot(data, label='data')
plt.legend()
modelp.name = 'Poisson'
modelnb.name = 'NegBin'

mods = OrderedDict()
mods[modelp]= tracep
mods[modelnb]= tracenb

comp = pm.stats.compare(mods)
comp
pm.plots.compareplot(comp)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1469248291571754, 'fraction_numerical': 0.037965072133637055, 'mean_word_length': 5.3647342995169085, 'pattern_counts': {'":': 0, '<': 9, '<?xml version=': 0, '>': 9, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '2', 'dir': 'jupyter-structured-clean-dedup', 'id': '604750', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 1147, 'n_tokens_neox': 1052, 'n_words': 203, 'path': 'ModelSelection/ModelComparison1.ipynb', 'repo_name': 'balarsen/pymc_learning', 'size': '2634'}