Python autograd.numpy.copy() Examples

def gbrbm_perturb(var_perturb_B, dx=50, dh=10):
    """
    Get a Gaussian-Bernoulli RBM problem where the first entry of the B matrix
    (the matrix linking the latent and the observation) is perturbed.

    - var_perturb_B: Gaussian noise variance for perturbing B.
    - dx: observed dimension
    - dh: latent dimension

    Return p (density), data source
    """
    with util.NumpySeedContext(seed=10):
        B = np.random.randint(0, 2, (dx, dh))*2 - 1.0
        b = np.random.randn(dx)
        c = np.random.randn(dh)
        p = density.GaussBernRBM(B, b, c)

        B_perturb = np.copy(B)
        if var_perturb_B > 1e-7:
            B_perturb[0, 0] = B_perturb[0, 0] + \
                np.random.randn(1)*np.sqrt(var_perturb_B)
        ds = data.DSGaussBernRBM(B_perturb, b, c, burnin=2000)

    return p, ds 

def _set_transmat_prior(self, transmat_prior_val):
        # new val needs be n_unique x n_unique
        # internally, n_components x n_components
        # _ntied_transmat_prior is
        # called to get n_components x n_components
        transmat_prior_new = np.zeros((self.n_components, self.n_components))
        if transmat_prior_val is not None:

            if transmat_prior_val.shape == (self.n_unique, self.n_unique):
                transmat_prior_new = \
                np.copy(self._ntied_transmat_prior(transmat_prior_val))

            else:
                raise ValueError("cannot match shape of transmat_prior")


        self.transmat_prior = transmat_prior_new 

def _set_transmat(self, transmat_val):
        if transmat_val is None:
            transmat = np.tile(1.0 / self.n_components,
                               (self.n_components, self.n_components))
        else:
            transmat_val[np.isnan(transmat_val)] = 0.0
            normalize(transmat_val, axis=1)

            if (np.asarray(transmat_val).shape == (self.n_components,
                                                   self.n_components)):
                transmat = np.copy(transmat_val)
            elif transmat_val.shape[0] == self.n_unique:
                transmat = self._ntied_transmat(transmat_val)
            else:
                raise ValueError("cannot match shape of transmat")

        if not np.all(np.allclose(np.sum(transmat, axis=1), 1.0)):
            raise ValueError('Rows of transmat must sum to 1.0')
        self._log_transmat = np.log(np.asarray(transmat).copy())
        underflow_idx = np.isnan(self._log_transmat)
        self._log_transmat[underflow_idx] = NEGINF 

def _set_startprob_prior(self, startprob_prior):
        if startprob_prior is None or startprob_prior == 1.0:
            startprob_prior = np.zeros(self.n_components)
        else:
            startprob_prior = np.asarray(startprob_prior, dtype=np.float)

            if len(startprob_prior) != self.n_components:
                if len(startprob_prior) == self.n_unique:
                    startprob_prior_split = np.copy(startprob_prior) / \
                        (1.0 + self.n_tied)
                    startprob_prior = np.zeros(self.n_components)
                    for u in range(self.n_unique):
                        for t in range(self.n_chain):
                            startprob_prior[u*(self.n_chain)+t] = \
                                startprob_prior_split[u].copy()
                else:
                    raise ValueError("cannot match shape of startprob")

        self.startprob_prior = np.asarray(startprob_prior).copy() 

def _set_startprob(self, startprob):
        if startprob is None:
            startprob = np.tile(1.0 / self.n_components, self.n_components)
        else:
            startprob = np.asarray(startprob, dtype=np.float)

            normalize(startprob)

            if len(startprob) != self.n_components:
                if len(startprob) == self.n_unique:
                    startprob_split = np.copy(startprob) / (1.0+self.n_tied)
                    startprob = np.zeros(self.n_components)
                    for u in range(self.n_unique):
                        for t in range(self.n_chain):
                            startprob[u*(self.n_chain)+t] = \
                                startprob_split[u].copy()
                else:
                    raise ValueError("cannot match shape of startprob")

        if not np.allclose(np.sum(startprob), 1.0):
            raise ValueError('startprob must sum to 1.0')

        self._log_startprob = np.log(np.asarray(startprob).copy()) 

def _set_precision(self, precision_val):
        # new val needs to have dimension (n_unique, n_features, n_features)
        # internally, n_components x 1
        precision_new = \
        np.zeros((self.n_components, self.n_features, self.n_features))
        if precision_val is not None:
            precision_val = \
            precision_val.reshape(self.n_unique, self.n_features, self.n_features)
            if precision_val.shape == \
            (self.n_unique, self.n_features, self.n_features):
                for u in range(self.n_unique):
                    for t in range(self.n_chain):
                        precision_new[u*(self.n_chain)+t] = precision_val[u].copy()
            else:
                raise ValueError("cannot match shape of precision")
        self._precision_ = precision_new 

def _set_alpha(self, alpha_val):
        # new val needs to have a 1st dim of length n_unique x n_lags
        # if shared_alpha is true, a shape of 1 x n_lags is possible, too
        # internally, n_components x n_lags
        alpha_new = np.zeros((self.n_components, self.n_lags))

        if alpha_val is not None:
            if alpha_val.ndim == 1:
                alpha_val = alpha_val.reshape(-1, 1)  # make sure 2nd dim exists

            if alpha_val.shape[1] != self.n_lags:
                raise ValueError("shape[1] does not match n_lags")

            if self.shared_alpha == False:
                # alpha is not shared
                if alpha_val.shape[0] != self.n_unique:
                    raise ValueError("shape[0] does not match n_unique")
                for u in range(self.n_unique):
                    for t in range(1+self.n_tied):
                        alpha_new[u*(1+self.n_tied)+t, :] = alpha_val[u, :].copy()
            else:
                # alpha is shared ...
                if alpha_val.shape[0] != self.n_unique and \
                  alpha_val.shape[0] != 1:
                    # ... the shape should either be 1 x L or U x L
                    raise ValueError("shape[0] is neither 1 nor does it match n_unique")
                if alpha_val.shape[0] == self.n_unique and \
                  not (alpha_val == alpha_val[0,:]).all():
                    # .. in case of U x L the rows need to be identical
                    raise ValueError("rows not identical (shared_alpha)")
                for u in range(self.n_unique):
                    for t in range(1+self.n_tied):
                        alpha_new[u*(1+self.n_tied)+t, :] = alpha_val[0, :].copy()

        self._alpha_ = alpha_new 

def _obj_grad(self, wrt, m, p, a, xn, xln, gn, entries='all', **kwargs):
        m = m.reshape(self.n_unique, self.n_features, 1)  # tm

        if wrt == 'm':
            wrt_num = 0
        elif wrt == 'p':
            wrt_num = 1
        elif wrt == 'a':
            wrt_num = 2
        else:
            raise ValueError('unknown parameter')
        res = grad(self._obj, wrt_num)(m, p, a, xn, xln, gn)

        if wrt == 'p' and self.n_features > 1:
            if entries == 'diag':
                res_new = \
                np.zeros((self.n_unique, self.n_features, self.n_features))
                for u in range(self.n_unique):
                    for f in range(self.n_features):
                        res_new[u,f,f] = res[u,f,f]
                res = np.copy(res_new)

            elif entries == 'offdiag':
                for u in range(self.n_unique):
                    for f in range(self.n_features):
                        res[u,f,f] = 0.

        res = np.array([res])
        return res 

def _set_mu_prior(self, mu_prior):
        if mu_prior is None:
            self._mu_prior_ = np.zeros((self.n_components, self.n_features))
        else:
            mu_prior = np.asarray(mu_prior)
            mu_prior = mu_prior.reshape(self.n_unique, self.n_features)
            if mu_prior.shape == (self.n_unique, self.n_features):
                for u in range(self.n_unique):
                    for t in range(self.n_chain):
                        self._mu_prior[u*(self.n_chain)+t] = mu_prior[u].copy()

            else:
                raise ValueError("cannot match shape of mu_prior") 

def _set_mu(self, mu_val):
        # new val needs to be of shape (n_uniqe, n_features)
        # internally, (n_components x n_features)
        mu_new = np.zeros((self.n_components, self.n_features))
        if mu_val is not None:
            mu_val = mu_val.reshape(self.n_unique, self.n_features)
            if mu_val.shape == (self.n_unique, self.n_features):
                for u in range(self.n_unique):
                    for t in range(self.n_chain):
                        mu_new[u*(self.n_chain)+t] = mu_val[u].copy()
            else:
                raise ValueError("cannot match shape of mu")
        self._mu_ = mu_new 

def _calc_pareto_front(self, n_pareto_points=100):
        x1 = anp.linspace(0, 5, n_pareto_points)
        x2 = anp.copy(x1)
        x2[x1 >= 3] = 3
        return anp.vstack((4 * anp.square(x1) + 4 * anp.square(x2), anp.square(x1 - 5) + anp.square(x2 - 5))).T 

def __uniform_reference_directions(self, ref_dirs, ref_dir, n_partitions, beta, depth):
        if depth == len(ref_dir) - 1:
            ref_dir[depth] = beta / (1.0 * n_partitions)
            ref_dirs.append(ref_dir[None, :])
        else:
            for i in range(beta + 1):
                ref_dir[depth] = 1.0 * i / (1.0 * n_partitions)
                self.__uniform_reference_directions(ref_dirs, anp.copy(ref_dir), n_partitions, beta - i,
                                                    depth + 1) 

def __add__(self, data2):
        """
        Merge the current Data with another one.
        Create a new Data and create a new copy for all internal variables.
        """
        copy = self.clone()
        copy2 = data2.clone()
        nX = np.vstack((copy.X, copy2.X))
        return Data(nX)

### end Data class 

def clone(self):
        """
        Return a new Data object with a separate copy of each internal 
        variable, and with the same content.
        """
        nX = np.copy(self.X)
        return Data(nX) 

def sample(self, n_samples=2000, observed_states=None, random_state=None):
        """Generate random samples from the self.

        Parameters
        ----------
        n : int
            Number of samples to generate.

        observed_states : array
            If provided, states are not sampled.

        random_state: RandomState or an int seed
            A random number generator instance. If None is given, the
            object's random_state is used

        Returns
        -------
        samples : array_like, length (``n_samples``, ``n_features``)
                  List of samples

        states : array_like, shape (``n_samples``)
                 List of hidden states (accounting for tied states by giving
                 them the same index)
        """

        if random_state is None:
            random_state = self.random_state
        random_state = check_random_state(random_state)

        samples = np.zeros((n_samples, self.n_features))
        states = np.zeros(n_samples)

        if observed_states is None:
            startprob_pdf = np.exp(np.copy(self._log_startprob))
            startdist = stats.rv_discrete(name='custm',
                                      values=(np.arange(startprob_pdf.shape[0]),
                                                        startprob_pdf),
                                      seed=random_state)
            states[0] = startdist.rvs(size=1)[0]

            transmat_pdf = np.exp(np.copy(self._log_transmat))
            transmat_cdf = np.cumsum(transmat_pdf, 1)

            nrand = random_state.rand(n_samples)
            for idx in range(1,n_samples):
                newstate = (transmat_cdf[int(states[idx-1])] > nrand[idx-1]).argmax()
                states[idx] = newstate

        else:
            states = observed_states

        mu = np.copy(self._mu_)
        precision = np.copy(self._precision_)
        for idx in range(n_samples):
            mean_ = mu[states[idx]]

            covar_ = np.linalg.inv(precision[states[idx]])
            samples[idx] = multivariate_normal.rvs(mean=mean_, cov=covar_,
                                                   random_state=random_state)
        states = self._process_sequence(states)
        return samples, states 

def _init_params(self, data, lengths=None, params='stmp'):
        X = data['obs']

        if 's' in params:
            self.startprob_.fill(1.0 / self.n_components)

        if 't' in params or 'm' in params or 'p' in params:

            kmmod = cluster.KMeans(n_clusters=self.n_unique,
                                   random_state=self.random_state).fit(X)
            kmeans = kmmod.cluster_centers_

        if 't' in params:
            # TODO: estimate transitions from data (!) / consider n_tied=1
            if self.n_tied == 0:
                transmat = np.ones([self.n_components, self.n_components])
                np.fill_diagonal(transmat, 10.0)
                self.transmat_ = transmat  # .90 for self-transition

            else:
                transmat = np.zeros((self.n_components, self.n_components))
                transmat[range(self.n_components),
                         range(self.n_components)] = 100.0  # diagonal
                transmat[range(self.n_components-1),
                         range(1, self.n_components)] = 1.0  # diagonal + 1
                transmat[[r * (self.n_chain) - 1
                          for r in range(1, self.n_unique+1)
                          for c in range(self.n_unique-1)],
                         [c * (self.n_chain)
                          for r in range(self.n_unique)
                          for c in range(self.n_unique) if c != r]] = 1.0

                self.transmat_ = np.copy(transmat)

        if 'm' in params:
            mu_init = np.zeros((self.n_unique, self.n_features))
            for u in range(self.n_unique):
                for f in range(self.n_features):
                    mu_init[u][f] = kmeans[u, f]

            self.mu_ = np.copy(mu_init)

        if 'p' in params:
            precision_init = np.zeros((self.n_unique, self.n_features, self.n_features))
            for u in range(self.n_unique):
                if self.n_features == 1:
                    precision_init[u] = np.linalg.inv(np.cov(X[kmmod.labels_ == u], bias = 1))
                else:
                    precision_init[u] = np.linalg.inv(np.cov(np.transpose(X[kmmod.labels_ == u])))

            self.precision_ = np.copy(precision_init) 

def match(self, model_frame, diff_kernels=None, convolution="fft"):
        """Match the frame of `Blend` to the frame of this observation.

        The method sets up the mappings in spectral and spatial coordinates,
        which includes a spatial selection, computing PSF difference kernels
        and filter transformations.

        Parameters
        ---------
        model_frame: a `scarlet.Frame` instance
            The frame of `Blend` to match
        diff_kernels: array
            The difference kernel for each band.
            If `diff_kernels` is `None` then they are
            calculated automatically.
        convolution: str
            The type of convolution to use.
            - `real`: Use a real space convolution and gradient
            - `fft`: Use a DFT to do the convolution and gradient

        Returns
        -------
        None
        """
        self.model_frame = model_frame
        if model_frame.channels is not None and self.frame.channels is not None:
            channel_origin = list(model_frame.channels).index(self.frame.channels[0])
            self.frame.origin = (channel_origin, *self.frame.origin[1:])

        slices = overlapped_slices(self.frame, model_frame)
        self.slices_for_images = slices[0] # Slice of images to match the model
        self.slices_for_model = slices[1] #  Slice of model that overlaps with the observation

        # check dtype consistency
        if self.frame.dtype != model_frame.dtype:
            self.frame.dtype = model_frame.dtype
            self.images = self.images.copy().astype(model_frame.dtype)
            if type(self.weights) is np.ndarray:
                self.weights = self.weights.copy().astype(model_frame.dtype)

        # constrcut diff kernels
        if diff_kernels is None:
            self._diff_kernels = None
            if self.frame.psf is not model_frame.psf:
                assert self.frame.psf is not None and model_frame.psf is not None
                psf = fft.Fourier(self.frame.psf.update_dtype(model_frame.dtype).image)
                model_psf = fft.Fourier(
                    model_frame.psf.update_dtype(model_frame.dtype).image
                )
                self._diff_kernels = fft.match_psfs(psf, model_psf)
        else:
            if not isinstance(diff_kernels, fft.Fourier):
                diff_kernels = fft.Fourier(diff_kernels)
            self._diff_kernels = diff_kernels

        # initialize the filter window
        assert convolution in ["real", "fft"], "`convolution` must be either 'real' or 'fft'"
        self.convolution = convolution
        return self