Start, end, and threshold that clips the value of a time series the most, given a limitted "capacity" and "rate"

    Assumes that signal can be linearly interpolated between points (trapezoidal integration)

    Arguments:
      ts (TimeSeries): Time series to attempt to clip to as low a max value as possible
      capacity (float): Total "funds" or "energy" available for clipping (integrated area under time series)
      method (str): scipy optimization algorithm name, one of:
        'L-BFGS-B': Byrd, 1995, "A Limited Memory Algorithm for Bound Constrained Optimization"
        'TNC':      Truncated Newton in C, or Newton Conjugate-Gradient, each variable may be constrained with upper and lower bounds
        'COBYLA':   Constrained Optimization by Linear Approximation. Fortran implementation.
        'SLSQP':    Kraft, 1988, Sequential Least Squares Programming or Quadratic Programming, infinite bounds converted to large floats

    TODO:
      Bisection search for the optimal threshold.

    Returns:
      2-tuple: Timestamp of the start and end of the period of the maximum clipped integrated increase

    >>> t = ['2014-12-09T00:00', '2014-12-09T00:15', '2014-12-09T00:30', '2014-12-09T00:45', '2014-12-09T01:00', '2014-12-09T01:15', '2014-12-09T01:30', '2014-12-09T01:45']
    >>> import pandas as pd
    >>> ts = pd.Series([217, 234, 235, 231, 219, 219, 231, 232], index=pd.to_datetime(t))  # doctest: +ELLIPSIS, +NORMALIZE_WHITESPACE
    >>> clipping_params(ts, capacity=60000)['threshold']  # doctest: +ELLIPSIS, +NORMALIZE_WHITESPACE
    218.13...
    >>> clipping_params(ts, capacity=30000)['threshold']  # doctest: +ELLIPSIS, +NORMALIZE_WHITESPACE
    224.15358...

def clipping_params(ts, capacity=100, rate_limit=float('inf'), method=None, max_attempts=100):
    """Start, end, and threshold that clips the value of a time series the most, given a limitted "capacity" and "rate"

    Assumes that signal can be linearly interpolated between points (trapezoidal integration)

    Arguments:
      ts (TimeSeries): Time series to attempt to clip to as low a max value as possible
      capacity (float): Total "funds" or "energy" available for clipping (integrated area under time series)
      method (str): scipy optimization algorithm name, one of:
        'L-BFGS-B': Byrd, 1995, "A Limited Memory Algorithm for Bound Constrained Optimization"
        'TNC':      Truncated Newton in C, or Newton Conjugate-Gradient, each variable may be constrained with upper and lower bounds
        'COBYLA':   Constrained Optimization by Linear Approximation. Fortran implementation.
        'SLSQP':    Kraft, 1988, Sequential Least Squares Programming or Quadratic Programming, infinite bounds converted to large floats

    TODO:
      Bisection search for the optimal threshold.

    Returns:
      2-tuple: Timestamp of the start and end of the period of the maximum clipped integrated increase

    >>> t = ['2014-12-09T00:00', '2014-12-09T00:15', '2014-12-09T00:30', '2014-12-09T00:45', '2014-12-09T01:00', '2014-12-09T01:15', '2014-12-09T01:30', '2014-12-09T01:45']
    >>> import pandas as pd
    >>> ts = pd.Series([217, 234, 235, 231, 219, 219, 231, 232], index=pd.to_datetime(t))  # doctest: +ELLIPSIS, +NORMALIZE_WHITESPACE
    >>> clipping_params(ts, capacity=60000)['threshold']  # doctest: +ELLIPSIS, +NORMALIZE_WHITESPACE
    218.13...
    >>> clipping_params(ts, capacity=30000)['threshold']  # doctest: +ELLIPSIS, +NORMALIZE_WHITESPACE
    224.15358...
    """
    VALID_METHODS = ['L-BFGS-B', 'TNC', 'SLSQP', 'COBYLA']
    # print('in clipping params for ts.index={0} and method={1}'.format(ts.index[0], method))
    ts.index = ts.index.astype(np.int64)
    costs = []

    def cost_fun(x, *args):
        thresh = x[0]
        ts, capacity, bounds = args
        integral = clipped_area(ts, thresh=thresh)
        terms = np.array([(10. * (integral - capacity) / capacity) ** 2,
                        2. / 0.1**((bounds[0] - thresh) * capacity / bounds[0]),
                        2. / 0.1**((thresh - bounds[1]) * capacity / bounds[1]),
                        1.2 ** (integral / capacity)])
        return sum(terms)

    bounds = (ts.min(), ts.max())
    done, attempts = 0, 0
    thresh0 = bounds[0] + 0.5 * (bounds[1] - bounds[0])
    if not method or not method in VALID_METHODS:
        while attempts < max_attempts and not done:
            for optimizer_method in VALID_METHODS:
                optimum = minimize(fun=cost_fun, x0=[thresh0], bounds=[bounds], args=(ts, capacity, bounds), method=optimizer_method)
                if optimum.success:
                    done = True
                    break
            if done:
                break
            attempts += 1
            thresh0 = bounds[0] + random.random() * (bounds[1] - bounds[0])
    else:
        optimum = minimize(fun=cost_fun, x0=[thresh0], bounds=[bounds], args=(ts, capacity, bounds), method=method)
    thresh = optimum.x[0]
    integral = clipped_area(ts, thresh=thresh)
    params = dict(optimum)
    params.update({'costs': costs, 'threshold': thresh, 'initial_guess': thresh0, 'attempts': attempts,
                   'integral': integral, 'method': method})
    return params