Writing a new distribution (2)

Questions
8

Does this work?

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
  ...

Alternatively you might manage with just modeling a uniform x and applying the power law function deterministically statistics - python scipy.stats.powerlaw negative exponent - Stack Overflow