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