I’m going to write a power-law distribution ( y = x^-2) in which x comes from a uniform distribution in [0,100] interval. According to this topic I could write the following code:
The way you wrote, the model is not making use of the model variable x anywhere. In the logp function, x is the value the sampler suggests for s1 at each iteration and it has nothing to do with the x variable in the model. You could have called it anything else (the standard is value). I am not sure how you wanted to incorporate x in your model, so I cannot give you a solution.
Can you show how you are generating your data and how the 2 unobserved variables x, s1 are used?
I have the feeling you don’t really want/ need x (s1 follows a power law distribution with your given alpha already, except it also accepts negative values as mentioned above) but I need you to clarify it.
To elaborate, I express my codes in terms of Numpy library. This is my target function: y = x^-2/3., and I want to draw samples from it, so we have:
np.random.power(a = 1/3)
note that in Numpy, the power law function is: y = ax^(a-1) (we can ignore the coefficient a!)
The above code implements what I want except that the default range for x is [0,1]. I want to extend it to [0,100].
x in my problem is a physical quantity which can vary between 0 ,100 and I have to consider this condition
def power(alpha):
def logp(value):
return tt.switch(
tt.and_(value > 0, value < 100),
tt.log(alpha) + (a-1)*tt.log(value),
-np.inf,
)
return logp
with pm.Model() as m:
y = pm.DensityDist(`y`, power(1/3)). # an interval transform here would be more efficient than the switch statement
...
Unfortunately, I couldn’t sample from your code!! It comes with big divergences (despite my target_accept which was 0.99!). Could you complete your code and sample from it?
Because applying a transformation to a random variable Y = f(X) does not transform the distribution in the same way (P(Y) \neq f(P(X))). In your case, as you’ve raised a uniform distribution to a power, you will obtain something like a Beta distribution (and not a power-law distribution).