def lp_model(trace, slice=slice(None), order=8, mode="f", mirror=None,
             method="svd", full=False):
    # check for bad arguments
    if mode not in ["f", "b"]:
        raise ValueError("mode must be 'f' or 'b'")
    if method not in ['svd', 'qr', 'cholesky', 'tls', 'hsvd']:
        raise ValueError("Invalid method")
    if trace.ndim != 1:
        raise ValueError("trace must be a 1D array")

    x = trace[slice]    # extract region to use for finding LP coefficients

    if mirror is not None:  # make mirror image if requested
        x = make_mirror(x, mirror)

    # calculate LP coefficient and factor to find poles
    if method in ['svd', 'qr', 'cholseky', 'tls']:
        D, d = make_Dd(x, order, mode)   # form the LP equation elements
        a = find_lpc(D, d, method)      # find LP coefficients
        poles = find_roots(a, mode)     # find roots

    elif method == "hsvd":
        poles = find_lproots_hsvd(x, M=order, K=order, mode=mode, zmethod='sm')
    else:
        raise ValueError("Invalid method")

    # reverse poles if we have backward poles
    if mode == "b":
        poles = [1. / pole for pole in poles]

    # determind the damping factor and frequencies from the roots
    damp = [root2damp(pole) for pole in poles]
    freq = [root2freq(pole) for pole in poles]

    if full is False:
        return damp, freq

    # perform Least Squares fitting to determind amplitudes and phases.

    # We need to find a least squares solutions to:
    # z_0*b_0^0+z_1*b_1^0+.... = x_0
    # z_0*b_0^1+z_1*b_1^1+.... = x_1
    # ...
    # Where b is the LP roots (poles), x is the signal, and a_0 are unknown
    #
    # To solve this recast into Bz = x and solve for a using scipy.lstsq

    # SPEED
    # B is a Vandermonde matrix, this characteristic can be used to more
    # efficiently solve this least squares problem.

    # build the B matrix (a Vandermonde matrix) and solve for the coefficients
    poles = np.array(poles)
    B = np.row_stack([poles ** (i) for i in range(len(x))])
    z, resid, rank, s = np.linalg.lstsq(B, np.array(x))

    # Now the z_n = amp_n*exp(phase_n*i), use this to determind the amplitudes
    # and phases
    amp = [cof2amp(cof) for cof in z]
    phase = [cof2phase(cof) for cof in z]

    return damp, freq, amp, phase


###############################################################
# functions to determine signal parameters from LP parameters #
###############################################################
Use Linear Prediction to model 1D NMR time domain data.
dataset/ETHPy150Open jjhelmus/nmrglue/nmrglue/process/proc_lp.py/lp_model