Can not correctly sample high dimensional RVs from customized likelihood

I tried to sample RVs from a customized likelihood whose parameters are known. In other word, I want to do something like

p0 = 0.3 #true bernoulli parameter
obs = stats.bernoulli.rvs(p0, size=1000) #obs comes from bernoulli distribution

but I need to modify two things: (1) use a customized likelihood to replace bernoulli distribution, and (2) my RV has two dimension, with the first dimension to be continuous and the second to be discrete.

I tried the following, but got something strange.

import pandas as pd #To work with dataset
import numpy as np #Math library
import seaborn as sns
import matplotlib.pyplot as plt #to plot some parameters
import pymc3 as pm
import theano.tensor as tt

	def my_density_1(theta,W):
		
	def logp_1(X,Y):
		
		def log_expo(lam,x):
			return( tt.sum(tt.log(lam) - lam * x) )

		def log_bernoulli(p,x):
			return( tt.sum(1.* ( tt.switch( x, tt.log(p), tt.log(1 - p) ))) )
		
		
		LL = np.array([
			log_expo(tetha[0],X), #first dimension, exponential RV
			log_bernoulli(tetha[1],Y), #first dimension, Bernoulli RV
			])
		
		return( np.dot(W,LL) ) # My customized likelihood: weighted sum of two components
	return(logp_1)

The above is my definition of cuctomized likelihood. The first component is a exponential distribution and the second component is a bernoulli distribution.
The RV is a two dimensional vector [X,Y]. My target is to get samples from this vector.
theta and W are parameters for the likelihood and they are known. Pleaase note I do NOT want to get posterior samples of these parameters. What I need is samples of random vector [X,Y].

Next I set pymc3 model using the following codes:

with pm.Model() as m:
	theta = np.array([
				3.384, #parametr for X
				0.397, #parametr for Y
				])
            W = [0.5,0.5] #weights of two components

	##set likelihood
	simulate = pm.DensityDist('simulate',my_density_1(theta,W))

When I run the above codes, I got the following error:

TypeError: logp_1() missing 1 required positional argument: 'Y'

What is wrong here?

Hello,

I encountered a similar error, did you figure out what might happen in your code? My code is something like below:

        with pm.Model() as m:
            sigma = pm.Uniform('sigma',lower=0,upper=1)
            nc = pm.HalfNormal('nc',sigma=sigma,observed=y)
            nc_hat = pm.Deterministic('nc_hat',pm.math.sum(nc)/len(y))
            psi = pm.Beta('alpha',alpha=2,beta=nc_hat*10)
            mu = pm.Gamma('beta',alpha=nc_hat*10,beta=1)
            c = pm.ZeroInflatedPoisson('c',psi=psi,mu=mu,observed=x)
        with m:
            step = pm.NUTS()
            trace = pm.sample(100,step=step,return_inferencedata=False,cores=1)

Error:

  File "/data/salomonis2/software/snaf/gtex.py", line 369, in tumor_specificity
    c = pm.ZeroInflatedPoisson('c',psi=psi,mu=mu,observed=x)
  File "/data/salomonis2/LabFiles/Frank-Li/refactor/neo_env/lib/python3.7/site-packages/pymc3/distributions/distribution.py", line 121, in __new__
    dist = cls.dist(*args, **kwargs)
  File "/data/salomonis2/LabFiles/Frank-Li/refactor/neo_env/lib/python3.7/site-packages/pymc3/distributions/distribution.py", line 130, in dist
    dist.__init__(*args, **kwargs)
TypeError: __init__() missing 1 required positional argument: 'theta'

Just curious, if you simply write

c = pm.ZeroInflatedPoisson('c',psi,mu,observed=x)

would something be changed?

1 Like

Mu is now called theta for consistency with the non zero inflated distribution (I think).

The original issue is completely different

1 Like

yes it works, thanks a lot

I see, thanks a lot for explaining it!