Changepoint detection with Beta MLE

beta_log_like.py

# WIP: This code does NOT work

import numpy as np
import matplotlib.pyplot as plt
from scipy.special import beta
from scipy.optimize import minimize


def neg_beta_log_likelihood(z, params):
    a0, b0, a1, b1, tau = params
    taui = int(np.round(tau))
    T = len(z)
    nll = -(
        (a0 - 1) * np.sum(np.log(z[:taui])) + (a1 - 1) * np.sum(np.log(z[taui:])) +
        + (b0 - 1) * np.sum(np.log(1 - z[:taui])) + (b1 - 1) * np.sum(np.log(1 - z[taui:])) 
        - (tau - 1) * np.log(beta(a0, b0)) - (T - tau) * np.log(beta(a1, b1)) 
    )
    return nll
     
def constraint(params):
    a0, b0, a1, b1, _ = params
    return - a0 / (a0 + b0) + a1 / (a1 + b1)

    
T = 10
z = np.hstack([
    np.random.beta(2, 10, size=T//2),
    np.random.beta(8, 10, size=T//2),   
])

res = minimize(
    fun = lambda params: neg_beta_log_likelihood(z, params),
    x0 = [2, 2, 2, 2, T//3],
    bounds = [(1, None)] * 4 + [(1, T-1)],
    constraints = [{'type': 'ineq', 'fun': constraint}],
    options={"maxiter":100000,  "ftol":1e-10,},
    method="SLSQP",
    tol=1e-10,
)
print(res)