<SYSTEM_TASK:>
Solve the following problem using Python, implementing the functions described below, one line at a time
<END_TASK>
<USER_TASK:>
Description:
 def solve_potts_autogamma(y, w, beta=None, **kw):
 """Solve Potts problem with automatically determined gamma. The optimal value is determined by minimizing the information measure:: f(gamma) = beta J(x(gamma)) + log sum(abs(x(gamma) - y)**p) where x(gamma) is the solution to the Potts problem for a fixed gamma. The minimization is only performed rather roughly. Parameters beta : float or 'bic' Penalty parameter. Default is 4*ln(n)/n, similar to Bayesian information criterion for gaussian model with unknown variance assuming 4 DOF per breakpoint. """

    n = len(y)

    if n == 0:
        return [], [], [], None

    mu_dist = get_mu_dist(y, w)
    mu, dist = mu_dist.mu, mu_dist.dist

    if beta is None:
        beta = 4 * math.log(n) / n

    gamma_0 = dist(0, n-1)

    if gamma_0 == 0:
        # Zero variance
        gamma_0 = 1.0

    best_r = [None]
    best_v = [None]
    best_d = [None]
    best_obj = [float('inf')]
    best_gamma = [None]

    def f(x):
        gamma = gamma_0 * math.exp(x)
        r, v, d = solve_potts_approx(y, w, gamma=gamma, mu_dist=mu_dist, **kw)

        # MLE fit noise correlation
        def sigma_star(rights, values, rho):
            """
            |E_0| + sum_{j>0} |E_j - rho E_{j-1}|
            """
            l = 1
            E_prev = y[0] - values[0]
            s = abs(E_prev)
            for r, v in zip(rights, values):
                for yv in y[l:r]:
                    E = yv - v
                    s += abs(E - rho*E_prev)
                    E_prev = E
                l = r
            return s

        rho_best = golden_search(lambda rho: sigma_star(r, v, rho), -1, 1,
                                 xatol=0.05, expand_bounds=True)

        # Measurement noise floor
        if len(v) > 2:
            absdiff = [abs(v[j+1] - v[j]) for j in range(len(v) - 1)]
            sigma_0 = 0.1 * min(absdiff)
        else:
            absv = [abs(z) for z in v]
            sigma_0 = 0.001 * min(absv)
        sigma_0 = max(1e-300, sigma_0)

        # Objective function
        s = sigma_star(r, v, rho_best)
        obj = beta*len(r) + math.log(sigma_0 + s)

        # Done
        if obj < best_obj[0]:
            best_r[0] = r
            best_v[0] = v
            best_d[0] = d
            best_gamma[0] = gamma
            best_obj[0] = obj
        return obj

    # Try to find best gamma (golden section search on log-scale); we
    # don't need an accurate value for it however
    a = math.log(0.1/n)
    b = 0.0
    golden_search(f, a, b, xatol=abs(a)*0.1, ftol=0, expand_bounds=True)
    return best_r[0], best_v[0], best_d[0], best_gamma[0]