Python numpy.vander() Examples

def __call__(self, xnew):
        saveshape = np.shape(xnew)
        xnew = np.ravel(xnew)
        res = np.empty_like(xnew)
        mask = (xnew >= self.a) & (xnew <= self.b)
        res[~mask] = self.fill
        xx = xnew.compress(mask)
        indxs = np.searchsorted(self.breaks, xx)-1
        indxs = indxs.clip(0, len(self.breaks))
        pp = self.coeffs
        diff = xx - self.breaks.take(indxs)
        V = np.vander(diff, N=self.K)
        # values = np.diag(dot(V,pp[:,indxs]))
        values = array([dot(V[k, :], pp[:, indxs[k]]) for k in xrange(len(xx))])
        res[mask] = values
        res.shape = saveshape
        return res 

def _ppoly_eval_2(coeffs, breaks, xnew, fill=np.nan):
    """Evaluate piecewise polynomial manually (another way)"""
    a = breaks[0]
    b = breaks[-1]
    K = coeffs.shape[0]

    saveshape = np.shape(xnew)
    xnew = np.ravel(xnew)
    res = np.empty_like(xnew)
    mask = (xnew >= a) & (xnew <= b)
    res[~mask] = fill
    xx = xnew.compress(mask)
    indxs = np.searchsorted(breaks, xx)-1
    indxs = indxs.clip(0, len(breaks))
    pp = coeffs
    diff = xx - breaks.take(indxs)
    V = np.vander(diff, N=K)
    values = np.array([np.dot(V[k, :], pp[:, indxs[k]]) for k in xrange(len(xx))])
    res[mask] = values
    res.shape = saveshape
    return res 

def linest(*args, **kwargs):  # Excel reference: https://support.office.com/en-us/article/LINEST-function-84d7d0d9-6e50-4101-977a-fa7abf772b6d

    Y = list(args[0].values())
    X = list(args[1].values())

    if len(args) == 3:
        const = args[2]
        if isinstance(const,str):
            const = (const.lower() == "true")
    else:
        const = True

    degree = kwargs.get('degree',1)

    # build the vandermonde matrix
    A = np.vander(X, degree+1)

    if not const:
        # force the intercept to zero
        A[:,-1] = np.zeros((1,len(X)))

    # perform the fit
    (coefs, residuals, rank, sing_vals) = np.linalg.lstsq(A, Y)

    return coefs 

def Nordsieck_RKn(self,t0,y,sw0):
        s=self.number_of_steps
        H=(s-1)*self.H
        co_nord=[N.array([1./2,1.]),N.array([2./5,3./5,1.])]
        l=size(y,0)
        y0=y[0,:]
        yf=self.f(t0,y0,sw0)
        
        if l==3:
            co=N.array([co_nord[0]])
            nord_n=N.vander(co_nord[0],self.number_of_steps+1)
            b=y[1:]-y0-co.T*yf
            nord=Sc.solve(nord_n[0:2,0:2],b)
        elif l==4:
            co=N.array([co_nord[1]])
            nord_n=N.vander(co_nord[1],self.number_of_steps+1)
            b=y[1:]-y0-H*co.T*yf
            nord=Sc.solve(nord_n[0:3,0:3],b)
        nord=N.vstack((y0,H*yf,nord[::-1]))       
        return nord 

def vander(x, n=None):
    """
    Masked values in the input array result in rows of zeros.

    """
    _vander = np.vander(x, n)
    m = getmask(x)
    if m is not nomask:
        _vander[m] = 0
    return _vander 

def vander(x, n=None):
    """
    Masked values in the input array result in rows of zeros.

    """
    _vander = np.vander(x, n)
    m = getmask(x)
    if m is not nomask:
        _vander[m] = 0
    return _vander 

def vander(x, n=None):
    """
    Masked values in the input array result in rows of zeros.

    """
    _vander = np.vander(x, n)
    m = getmask(x)
    if m is not nomask:
        _vander[m] = 0
    return _vander 

def vander(x, n=None):
    """
    Masked values in the input array result in rows of zeros.

    """
    _vander = np.vander(x, n)
    m = getmask(x)
    if m is not nomask:
        _vander[m] = 0
    return _vander 

def vander(x, n=None):
    """
    Masked values in the input array result in rows of zeros.

    """
    _vander = np.vander(x, n)
    m = getmask(x)
    if m is not nomask:
        _vander[m] = 0
    return _vander 

def detrend(x, order=1, axis=0):
    """
    Detrend an array with a trend of given order along axis 0 or 1

    Parameters
    ----------
    x : array_like, 1d or 2d
        data, if 2d, then each row or column is independently detrended with the
        same trendorder, but independent trend estimates
    order : int
        specifies the polynomial order of the trend, zero is constant, one is
        linear trend, two is quadratic trend
    axis : int
        axis can be either 0, observations by rows,
        or 1, observations by columns

    Returns
    -------
    detrended data series : ndarray
        The detrended series is the residual of the linear regression of the
        data on the trend of given order.
    """
    if x.ndim == 2 and int(axis) == 1:
        x = x.T
    elif x.ndim > 2:
        raise NotImplementedError('x.ndim > 2 is not implemented until it is needed')

    nobs = x.shape[0]
    if order == 0:
        # Special case demean
        resid = x - x.mean(axis=0)
    else:
        trends = np.vander(np.arange(float(nobs)), N=order + 1)
        beta = np.linalg.pinv(trends).dot(x)
        resid = x - np.dot(trends, beta)

    if x.ndim == 2 and int(axis) == 1:
        resid = resid.T

    return resid 

def reset_ramsey(res, degree=5):
    '''Ramsey's RESET specification test for linear models

    This is a general specification test, for additional non-linear effects
    in a model.


    Notes
    -----
    The test fits an auxiliary OLS regression where the design matrix, exog,
    is augmented by powers 2 to degree of the fitted values. Then it performs
    an F-test whether these additional terms are significant.

    If the p-value of the f-test is below a threshold, e.g. 0.1, then this
    indicates that there might be additional non-linear effects in the model
    and that the linear model is mis-specified.


    References
    ----------
    http://en.wikipedia.org/wiki/Ramsey_RESET_test

    '''
    order = degree + 1
    k_vars = res.model.exog.shape[1]
    #vander without constant and x:
    y_fitted_vander = np.vander(res.fittedvalues, order)[:, :-2] #drop constant
    exog = np.column_stack((res.model.exog, y_fitted_vander))
    res_aux = OLS(res.model.endog, exog).fit()
    #r_matrix = np.eye(degree, exog.shape[1], k_vars)
    r_matrix = np.eye(degree-1, exog.shape[1], k_vars)
    #df1 = degree - 1
    #df2 = exog.shape[0] - degree - res.df_model  (without constant)
    return res_aux.f_test(r_matrix) #, r_matrix, res_aux 

def vander(x, n=None):
    """
    Masked values in the input array result in rows of zeros.
    """
    _vander = np.vander(x, n)
    m = getmask(x)
    if m is not nomask:
        _vander[m] = 0
    return _vander 

def test_graph_laplacian():
    mats = ('np.arange(10) * np.arange(10)[:, np.newaxis]',
            'np.ones((7, 7))',
            'np.eye(19)',
            'sparse.diags([1, 1], [-1, 1], shape=(4,4))',
            'sparse.diags([1, 1], [-1, 1], shape=(4,4)).todense()',
            'np.asarray(sparse.diags([1, 1], [-1, 1], shape=(4,4)).todense())',
            'np.vander(np.arange(4)) + np.vander(np.arange(4)).T',
            )

    for mat_str in mats:
        for normed in (True, False):
            yield _check_graph_laplacian, mat_str, normed 

def test_gmres_basic():
    A = np.vander(np.arange(10) + 1)[:, ::-1]
    b = np.zeros(10)
    b[0] = 1
    x = np.linalg.solve(A, b)

    x_gm, err = gmres(A, b, restart=5, maxiter=1)

    assert_allclose(x_gm[0], 0.359, rtol=1e-2) 

def vander(x, n=None):
    """
    Masked values in the input array result in rows of zeros.

    """
    _vander = np.vander(x, n)
    m = getmask(x)
    if m is not nomask:
        _vander[m] = 0
    return _vander 

def testVander(self, shape, dtype, n, increasing, rng_factory):
    rng = rng_factory()
    def onp_fun(arg):
      arg = arg.astype(onp.float32) if dtype == lnp.bfloat16 else arg
      return onp.vander(arg, N=n, increasing=increasing)
    lnp_fun = lambda arg: lnp.vander(arg, N=n, increasing=increasing)
    args_maker = lambda: [rng([shape], dtype)]
    # np.vander seems to return float64 for all floating types. We could obey
    # those semantics, but they seem like a bug.
    self._CheckAgainstNumpy(onp_fun, lnp_fun, args_maker, check_dtypes=False,
                            tol={onp.float32: 1e-3})
    self._CompileAndCheck(
        lnp_fun, args_maker, check_dtypes=False, check_incomplete_shape=True) 

def vander(x, N=None, increasing=False):  # pylint: disable=missing-docstring,invalid-name
  x = asarray(x).data

  x_shape = tf.shape(x)
  N = N or x_shape[0]

  N_temp = utils.get_static_value(N)  # pylint: disable=invalid-name
  if N_temp is not None:
    N = N_temp
    if N < 0:
      raise ValueError('N must be nonnegative')
  else:
    tf.debugging.Assert(N >= 0, [N])

  rank = tf.rank(x)
  rank_temp = utils.get_static_value(rank)
  if rank_temp is not None:
    rank = rank_temp
    if rank != 1:
      raise ValueError('x must be a one-dimensional array')
  else:
    tf.debugging.Assert(rank == 1, [rank])

  if increasing:
    start = 0
    limit = N
    delta = 1
  else:
    start = N - 1
    limit = -1
    delta = -1

  x = tf.expand_dims(x, -1)
  return utils.tensor_to_ndarray(
      tf.math.pow(x, tf.cast(tf.range(start, limit, delta), dtype=x.dtype))) 

def vander(x, n=None):
    """
    Masked values in the input array result in rows of zeros.

    """
    _vander = np.vander(x, n)
    m = getmask(x)
    if m is not nomask:
        _vander[m] = 0
    return _vander 

def test_symmetric_graph_laplacian():
    symmetric_mats = ('np.arange(10) * np.arange(10)[:, np.newaxis]',
            'np.ones((7, 7))',
            'np.eye(19)',
            'sparse.diags([1, 1], [-1, 1], shape=(4,4))',
            'sparse.diags([1, 1], [-1, 1], shape=(4,4)).todense()',
            'np.asarray(sparse.diags([1, 1], [-1, 1], shape=(4,4)).todense())',
            'np.vander(np.arange(4)) + np.vander(np.arange(4)).T')
    for mat_str in symmetric_mats:
        for normed in True, False:
            _check_symmetric_graph_laplacian(mat_str, normed) 

def test_gmres_basic():
    A = np.vander(np.arange(10) + 1)[:, ::-1]
    b = np.zeros(10)
    b[0] = 1
    x = np.linalg.solve(A, b)

    with suppress_warnings() as sup:
        sup.filter(DeprecationWarning, ".*called without specifying.*")
        x_gm, err = gmres(A, b, restart=5, maxiter=1)

    assert_allclose(x_gm[0], 0.359, rtol=1e-2) 

def vander(x, n=None):
    """
    Masked values in the input array result in rows of zeros.

    """
    _vander = np.vander(x, n)
    m = getmask(x)
    if m is not nomask:
        _vander[m] = 0
    return _vander 

def vander(x, n=None):
    """
    Masked values in the input array result in rows of zeros.
    """
    _vander = np.vander(x, n)
    m = getmask(x)
    if m is not nomask:
        _vander[m] = 0
    return _vander 

def vander(x, n=None):
    """
    Masked values in the input array result in rows of zeros.

    """
    _vander = np.vander(x, n)
    m = getmask(x)
    if m is not nomask:
        _vander[m] = 0
    return _vander 

def vander(x, n=None):
    """
    Masked values in the input array result in rows of zeros.

    """
    _vander = np.vander(x, n)
    m = getmask(x)
    if m is not nomask:
        _vander[m] = 0
    return _vander 

def albrecht_5():
    # The values are solutions of
    # 6317094x^3 - 10022245*x^2 + 4149900*x - 336375 = 0
    sigma2 = roots([6317094, -10022245, 4149900, -336375])
    A = numpy.vander(sigma2, increasing=True).T
    b = numpy.array([frac(168899, 1350000), frac(7661, 180000), frac(71, 3000)])
    B = linear_solve(A, b)

    sqrt19 = sqrt(19)

    # ERR Stroud incorrectly lists sqrt(10) for s1.
    s1, s2 = sqrt((125 - pm_ * 10 * sqrt19) / 366)

    # ERR Stroud incorrectly lists 749489_3_.0 instead of 749489_2_.0
    C1, C2 = (7494892 + pm_ * 1053263 * sqrt19) / 205200000
    D = frac(81, 3125)

    u = sqrt(frac(5, 6)) * cos(pi / 8)
    v = sqrt(frac(5, 6)) * sin(pi / 8)

    data = [
        (B[0], fsd(2, (sqrt(sigma2[0]), 1))),
        (B[1], fsd(2, (sqrt(sigma2[1]), 1))),
        (B[2], fsd(2, (sqrt(sigma2[2]), 1))),
        (C1, pm([s1, s1])),
        (C2, pm([s2, s2])),
        (D, fsd(2, (u, 1), (v, 1))),
    ]

    points, weights = untangle(data)
    return S2Scheme("Albrecht 5", weights, points, 11, _source) 

def albrecht_6():
    # The values are solutions of
    # 11025*x^3 - 19020*x^2 + 9370*x - 1212 = 0
    sigma2 = roots([11025, -19020, 9370, -1212])
    A = numpy.vander(sigma2, increasing=True).T
    b = numpy.array([frac(1432433, 18849024), frac(1075, 31104), frac(521, 25920)])
    B = linear_solve(A, b)

    B0 = frac(2615, 43632)
    C = frac(16807, 933120)

    alpha = 2 * numpy.arange(10) * pi / 10
    rs = numpy.array([cos(alpha), sin(alpha)]).T

    alpha = (2 * numpy.arange(10) + 1) * pi / 10
    uv = numpy.array([cos(alpha), sin(alpha)]).T

    data = [
        (B0, z(2)),
        (B[0], sqrt(sigma2[0]) * rs),
        (B[1], sqrt(sigma2[1]) * rs),
        (B[2], sqrt(sigma2[2]) * rs),
        (C, sqrt(frac(6, 7)) * uv),
    ]

    points, weights = untangle(data)
    return S2Scheme("Albrecht 6", weights, points, 13, _source) 

def albrecht_7():
    alpha = 2 * numpy.arange(8) * pi / 8
    s = numpy.array([cos(alpha), sin(alpha)]).T

    alpha = (2 * numpy.arange(8) + 1) * pi / 8
    t = numpy.array([cos(alpha), sin(alpha)]).T

    sqrt21 = sqrt(21)
    wt1, wt2 = (4998 + pm_ * 343 * sqrt21) / 253125
    tau1, tau2 = sqrt((21 - pm_ * sqrt21) / 28)

    # The values are solutions of
    # 4960228*x^4 - 10267740*x^3 + 6746490*x^2 - 1476540*x + 70425 = 0
    sigma2 = roots([4960228, -10267740, 6746490, -1476540, 70425])
    A = numpy.vander(sigma2, increasing=True).T
    b = numpy.array(
        [frac(57719, 675000), frac(9427, 270000), frac(193, 9000), frac(113, 7200)]
    )
    ws = linear_solve(A, b)

    data = [
        (ws[0], sqrt(sigma2[0]) * s),
        (ws[1], sqrt(sigma2[1]) * s),
        (ws[2], sqrt(sigma2[2]) * s),
        (ws[3], sqrt(sigma2[3]) * s),
        (wt1, tau1 * t),
        (wt2, tau2 * t),
    ]

    points, weights = untangle(data)
    return S2Scheme("Albrecht 7", weights, points, 15, _source) 

def albrecht_8():
    alpha = 2 * numpy.arange(10) * pi / 10
    s = numpy.array([cos(alpha), sin(alpha)]).T

    alpha = (2 * numpy.arange(10) + 1) * pi / 10
    t = numpy.array([cos(alpha), sin(alpha)]).T

    m0 = frac(496439663, 13349499975)

    sqrt7 = sqrt(7)
    wt1, wt2 = (125504 + pm_ * 16054 * sqrt7) / 8751645
    tau1, tau2 = sqrt((14 - pm_ * sqrt7) / 18)

    # The values are solutions of
    # 160901628*x^4 - 364759920*x^3 + 274856190*x^2 - 76570340*x
    # + 6054195 = 0
    sigma2 = roots([160901628, -364759920, 274856190, -76570340, 6054195])
    A = numpy.vander(sigma2, increasing=True).T
    b = numpy.array(
        [
            frac(121827491812, 1802182496625),
            frac(48541, 1666980),
            frac(977, 55566),
            frac(671, 52920),
        ]
    )
    ws = linear_solve(A, b)

    data = [
        (m0, z(2)),
        (ws[0], sqrt(sigma2[0]) * s),
        (ws[1], sqrt(sigma2[1]) * s),
        (ws[2], sqrt(sigma2[2]) * s),
        (ws[3], sqrt(sigma2[3]) * s),
        (wt1, tau1 * t),
        (wt2, tau2 * t),
    ]

    points, weights = untangle(data)
    return S2Scheme("Albrecht 8", weights, points, 17, _source) 

def vander(x, n=None):
    """
    Masked values in the input array result in rows of zeros.
    """
    _vander = np.vander(x, n)
    m = getmask(x)
    if m is not nomask:
        _vander[m] = 0
    return _vander 

def _setup(self, config):
        torch.manual_seed(config['seed'])
        self.model = ButterflyProduct(size=config['size'],
                                      complex=False,
                                      fixed_order=config['fixed_order'],
                                      softmax_fn=config['softmax_fn'])
        if (not config['fixed_order']) and config['softmax_fn'] == 'softmax':
            self.semantic_loss_weight = config['semantic_loss_weight']
        self.optimizer = optim.Adam(self.model.parameters(), lr=config['lr'])
        self.n_steps_per_epoch = config['n_steps_per_epoch']
        size = config['size']
        # Need to transpose as dct acts on rows of matrix np.eye, not columns
        n = size
        np.random.seed(0)
        x = np.random.randn(n)
        V = np.vander(x, increasing=True)
        self.target_matrix = torch.tensor(V, dtype=torch.float)
        arange_ = np.arange(size)
        dct_perm = np.concatenate((arange_[::2], arange_[::-2]))
        br_perm = bitreversal_permutation(size)
        assert config['perm'] in ['id', 'br', 'dct']
        if config['perm'] == 'id':
            self.perm = torch.arange(size)
        elif config['perm'] == 'br':
            self.perm = br_perm
        elif config['perm'] == 'dct':
            self.perm = torch.arange(size)[dct_perm][br_perm]
        else:
            assert False, 'Wrong perm in config'