How to implement Beta Function effectively in Pymc3

Thank you for your help @ricardoV94! I switched to log-beta function and then converted the list into a tensor using the correct call and that fixed everything. My final code was:

from pymc3.distributions.dist_math import betaln

with pm.Model() as sbg_model:
    alpha = pm.Uniform('alpha', 0.0001, 10) # MLE est was 0.668
    beta = pm.Uniform('beta', 0.0001, 10) # MLE est was 3.806
    
    def logp(survived, churned, alpha = alpha, beta = beta):   
        # Calculate the final surivival probability (eq. 6)
        ln_surv_prob = np.log(beta_func(alpha, beta + n_obs)/beta_func(alpha, beta))
        
        # Find the probability of churn for all values prior
        ln_prob_vec = []
        for i in range(1, n_obs + 1):
            ln_prob_vec.append(betaln(alpha + 1, beta + i - 1) - betaln(alpha, beta))
        ln_prob_vec = tt.as_tensor_variable(ln_prob_vec, 'ln_prob_vec', ndim = 1)

        return pm.math.dot(churned, ln_prob_vec) + survived[-1] * ln_surv_prob   
    
    likelihood = pm.DensityDist('likelihood', logp, observed = {'survived': survived, 'churned': churned})
    sbg_trace = pm.sample(200) # just a sanity check to see if it even runs...
2 Likes