Fitting beta-binomial with different sample size (different n)

Questions v5
1

Hi everyone,

I have learned that the pm.Binomial function can take a tensor as an argument for n, and the observed data can also be a tensor. Below is an example code and a graph of my model. After the inference was completed, I generated posterior predictive samples, resulting in an x-array with shape (chain, draw, len(ad_ar)). Each draw was obtained from a posterior sample, and an array with shape (len(ad_ar),)… My question is whether the first element of this array was sampled from a Binomial distribution with parameters Binom(dp_ar[0], p). In other words, I want to know whether the posterior predictive samples were generated from multiple binomial distributions with different values for n (in this case, n was kept in the array ‘dp_ar’) when calculating the log-likelihood or posterior predictive samples.

ad_ar = np.array([229., 226., 403., 105., 295., 147., 486., 231., 330., 208., 246.,
       302., 104., 447., 161., 181., 210., 224., 288., 222., 172., 328.,
       342., 474., 300., 411., 337., 341., 250., 454., 506., 495., 156.,
       234., 261., 435., 252., 378., 309., 376., 374., 370., 277., 375.,
       327., 356., 454., 216., 286., 147., 177., 233., 152., 426., 399.,
       338., 280., 166., 113., 385., 271., 360., 283., 116., 243., 164.,
       311., 148., 208., 346., 224., 402., 436., 408., 199., 414., 323.,
       172., 234., 116., 183., 125., 385., 452., 218., 290., 173., 380.,
       459., 387., 289., 156., 518., 248., 185., 373., 301., 296., 356.,
       245., 165., 310., 215., 477., 517., 334., 636., 177., 241., 472.,
       201., 193., 332., 535., 306., 229., 433., 282., 353., 179., 322.,
       436., 292., 185., 485., 362., 228., 208., 310., 293., 358., 189.,
       105., 383., 195., 457., 433., 355., 221., 144., 197., 287., 420.,
       362., 327., 195., 430., 284., 174., 319., 411., 289., 295., 448.,
       183., 393., 458., 132., 183., 196., 508., 124., 343., 221., 365.,
       195., 100., 338., 207., 286., 282., 214., 168., 411., 169., 348.,
       150., 337., 354., 404., 398., 353., 127., 320., 395., 341., 271.,
       258., 234., 335., 453., 364., 202.])

dp_ar = np.array([433., 396., 403., 195., 295., 267., 486., 398., 330., 366., 452.,
       302., 191., 447., 293., 319., 378., 414., 288., 382., 370., 330.,
       342., 474., 300., 412., 337., 341., 250., 454., 506., 495., 263.,
       405., 261., 435., 252., 378., 309., 376., 374., 371., 278., 375.,
       328., 356., 455., 414., 289., 235., 306., 233., 261., 426., 399.,
       338., 280., 304., 224., 385., 271., 360., 283., 237., 420., 306.,
       311., 272., 369., 347., 387., 402., 437., 408., 366., 414., 323.,
       330., 430., 221., 358., 238., 385., 452., 218., 290., 288., 381.,
       459., 387., 290., 311., 518., 248., 357., 373., 301., 297., 356.,
       246., 300., 310., 370., 478., 517., 334., 637., 340., 241., 472.,
       375., 393., 332., 535., 306., 408., 433., 282., 353., 302., 322.,
       437., 293., 329., 485., 362., 421., 387., 311., 294., 358., 353.,
       232., 383., 373., 457., 433., 355., 385., 290., 340., 287., 420.,
       362., 327., 363., 431., 284., 296., 319., 411., 289., 296., 449.,
       184., 393., 458., 216., 320., 394., 508., 223., 343., 425., 365.,
       393., 217., 338., 371., 286., 282., 413., 296., 411., 304., 349.,
       297., 339., 354., 404., 399., 353., 230., 322., 395., 345., 271.,
       259., 235., 335., 454., 365., 370.])


beta_binom = pm.Model()
with beta_binom:
    p = pm.Beta("p", alpha=1, beta=1)
    ad = pm.Binomial("ad", n= dp_ar , p=p, observed=ad_ar)
    trace_beta_binom = pm.sample(2000,chains=4)
    pm.compute_log_likelihood(trace_beta_binom)
    pm.sample_posterior_predictive(trace_beta_binom,extend_inferencedata=True,random_seed=2023)

image

image

2

Yes, you get one random sample/ logp term for each entry in n.

You can read a bit about parameter dimensions in Distribution Dimensionality — PyMC 5.5.0 documentation

Everything about samples applies to probabilities, as you could always set samples as observed and the same shape rules applies.

1 Like
3

Thank you very much for your response. I have another question I hope you can help me with: from comparing the posterior predictive distribution with the observed data, it seems that the model is not fitting very well, with the observed data distribution being noticeably more discrete. Initially, I thought that this was due to generating posterior predictive samples that did not match the corresponding experiment’s n (the parameter of the binomial distribution) size, and instead used the same n parameter for all experiments, leading to this situation. However, it now appears that this is not the reason. So, what could be a possible reason for my model not being able to fit the data, even when I switch to the negative binomial distribution, with similar results? Is it possible that the sample experiments from these different groups are not generated from the same p (Binom(n,p)) samples?