Calculate the covariance matrix.

        Parameters
        ----------
        profile : site.Profile
            Input site profile

        Yields
        ------
        covar : `class`:numpy.array
            Covariance matrix
770538
def _calc_covar_matrix(self, profile):
        
        corr = self._calc_corr(profile)
        std = self._calc_ln_std(profile)
        # Modify the standard deviation by the truncated norm scale
        std *= randnorm.scale

        var = std ** 2
        covar = corr * std[:-1] * std[1:]

        # Main diagonal is the variance
        mat = diags([covar, var, covar], [-1, 0, 1]).toarray()

        return mat
1120495
def distribute_covar_matrix_to_match_covariance_type(
        tied_cv, covariance_type, n_components):
    """"""
    if covariance_type == 'spherical':
        cv = np.tile(tied_cv.mean() * np.ones(tied_cv.shape[1]),
                     (n_components, 1))
    elif covariance_type == 'tied':
        cv = tied_cv
    elif covariance_type == 'diag':
        cv = np.tile(np.diag(tied_cv), (n_components, 1))
    elif covariance_type == 'full':
        cv = np.tile(tied_cv, (n_components, 1, 1))
    else:
        raise ValueError("covariance_type must be one of " +
                         "'spherical', 'tied', 'diag', 'full'")
    return cv
188936
def fitlin(imgarr,refarr):
    
    # Initialize variables
    _mat = np.zeros((3,3),dtype=np.float64)
    _xorg = imgarr[0][0]
    _yorg = imgarr[0][1]
    _xoorg = refarr[0][0]
    _yoorg = refarr[0][1]
    _sigxox = 0.
    _sigxoy = 0.
    _sigxo = 0.
    _sigyox = 0.
    _sigyoy = 0.
    _sigyo = 0.

    _npos = len(imgarr)
    # Populate matrices
    for i in range(_npos):
        _mat[0][0] += np.power((imgarr[i][0] - _xorg),2)
        _mat[0][1] += (imgarr[i][0] - _xorg) * (imgarr[i][1] - _yorg)
        _mat[0][2] += (imgarr[i][0] - _xorg)
        _mat[1][1] += np.power((imgarr[i][1] - _yorg),2)
        _mat[1][2] += imgarr[i][1] - _yorg

        _sigxox += (refarr[i][0] - _xoorg)*(imgarr[i][0] - _xorg)
        _sigxoy += (refarr[i][0] - _xoorg)*(imgarr[i][1] - _yorg)
        _sigxo += refarr[i][0] - _xoorg
        _sigyox += (refarr[i][1] - _yoorg)*(imgarr[i][0] -_xorg)
        _sigyoy += (refarr[i][1] - _yoorg)*(imgarr[i][1] - _yorg)
        _sigyo += refarr[i][1] - _yoorg

    _mat[2][2] = _npos
    _mat[1][0] = _mat[0][1]
    _mat[2][0] = _mat[0][2]
    _mat[2][1] = _mat[1][2]

    # Now invert this matrix
    _mat = linalg.inv(_mat)

    _a  = _sigxox*_mat[0][0]+_sigxoy*_mat[0][1]+_sigxo*_mat[0][2]
    _b  = -1*(_sigxox*_mat[1][0]+_sigxoy*_mat[1][1]+_sigxo*_mat[1][2])
    #_x0 = _sigxox*_mat[2][0]+_sigxoy*_mat[2][1]+_sigxo*_mat[2][2]

    _c  = _sigyox*_mat[1][0]+_sigyoy*_mat[1][1]+_sigyo*_mat[1][2]
    _d  = _sigyox*_mat[0][0]+_sigyoy*_mat[0][1]+_sigyo*_mat[0][2]
    #_y0 = _sigyox*_mat[2][0]+_sigyoy*_mat[2][1]+_sigyo*_mat[2][2]

    _xt = _xoorg - _a*_xorg+_b*_yorg
    _yt = _yoorg - _d*_xorg-_c*_yorg

    return [_a,_b,_xt],[_c,_d,_yt]
630441
def generate_covariance(ts, sigma, tau):

    
    """

    ndim = ts.shape[0]

    tis = ts[:, np.newaxis]
    tjs = ts[np.newaxis, :]

    return sigma*sigma*np.exp(-np.square(tis-tjs)/(2.0*tau*tau))
188201
def covariance(self):
        

        mu = self.moments_central
        if mu[0, 0] != 0:
            m = mu / mu[0, 0]
            covariance = self._check_covariance(
                np.array([[m[0, 2], m[1, 1]], [m[1, 1], m[2, 0]]]))
            return covariance * u.pix**2
        else:
            return np.empty((2, 2)) * np.nan * u.pix**2
941332
def readCovarianceMatrixFile(cfile,readCov=True,readEig=True):
    "
    covFile = cfile+'.cov'
    evalFile = cfile+'.cov.eval'
    evecFile = cfile+'.cov.evec'

    RV = {}
    if readCov:
        assert os.path.exists(covFile), '%s is missing.'%covFile
        RV['K'] = SP.loadtxt(covFile)
    if readEig:
        assert os.path.exists(evalFile), '%s is missing.'%evalFile
        assert os.path.exists(evecFile), '%s is missing.'%evecFile
        RV['eval'] = SP.loadtxt(evalFile)
        RV['evec'] = SP.loadtxt(evecFile)

    return RV
1346027
protected function calculateClassCov()
    {
        // s is an n (number of classes) by m (number of column) matrix
        $s = array_fill(0, count($this->overallMean), array_fill(0, count($this->overallMean), 0));
        $sB = new Matrix($s, false);

        foreach ($this->means as $index => $classMeans) {
            $row = $this->calculateVar($classMeans, $this->overallMean);
            $N = $this->counts[$index];
            $sB = $sB->add($row->multiplyByScalar($N));
        }

        return $sB;
    }
1283107
def execute(self, **minimize_options):
        
        minimizer_ans = self.minimizer.execute(**minimize_options)
        try: # to build covariance matrix
            cov_matrix = minimizer_ans.covariance_matrix
        except AttributeError:
            cov_matrix = self.covariance_matrix(dict(zip(self.model.params, minimizer_ans._popt)))
        else:
            if cov_matrix is None:
                cov_matrix = self.covariance_matrix(dict(zip(self.model.params, minimizer_ans._popt)))
        finally:
            minimizer_ans.covariance_matrix = cov_matrix
        # Overwrite the DummyModel with the current model
        minimizer_ans.model = self.model
        minimizer_ans.gof_qualifiers['r_squared'] = r_squared(self.model, minimizer_ans, self.data)
        return minimizer_ans
1319090
def noise_covariance(fit, dof=2, **kw):
    
    ev = fit.eigenvalues

    measurement_noise = ev[-1]/(fit.n-dof)
    return 4*ev*measurement_noise
1157845
def is_trajectory_profile_orthogonal(nc, variable):
    '''
    
    '''
    # x(i, o), y(i, o), z(z), t(i, o)
    # X(i, o, z)
    dims = nc.variables[variable].dimensions
    cmatrix = coordinate_dimension_matrix(nc)

    for req in ('x', 'y', 'z', 't'):
        if req not in cmatrix:
            return False

    if len(cmatrix['x']) != 2:
        return False
    if cmatrix['x'] != cmatrix['y']:
        return False
    if cmatrix['x'] != cmatrix['t']:
        return False

    i, o = cmatrix['x']

    z = get_z_variable(nc)
    if cmatrix['z'] != (z,):
        return False

    if dims == (i, o, z):
        return True
    return False
731013
def getCovariance(self,normalize=True,i0=None,i1=None,pos0=None,pos1=None,chrom=None,center=True,unit=True,pos_cum0=None,pos_cum1=None,blocksize=None,X=None,**kw_args):
        """"""
        if X is not None:
            K=X.dot(X.T)
            Nsnp=X.shape[1]
        else:
            if (i0 is None) and (i1 is None) and ((pos0 is not None) & (pos1 is not None) & (chrom is not None)) or ((pos_cum0 is not None) & (pos_cum1 is not None)):
                i0,i1=self.getGenoIndex(pos0=pos0,pos1=pos1,chrom=chrom,pos_cum0=pos_cum0,pos_cum1=pose_cum1)

            [N,M]=self.genoM.shape
            if blocksize is None:
                blocksize=M
            if i0 is None:
                i0=0
            if i1 is None:
                i1=M
            nread = i0
            K=None
            Nsnp=i1-i0
            while nread<i1:
                thisblock=min(blocksize,i1-nread)
                X=self.getGenotypes(i0=nread,i1=(nread+thisblock),center=center,unit=unit,**kw_args)    
                if K is None:
                    K=X.dot(X.T)
                else:
                    K+=X.dot(X.T)
                nread+=thisblock
        if normalize:
            K/=(K.diagonal().mean())
        else:#divide by number of SNPs in K
            K/=Nsnp
        return K
734187
private function computeCovariance()
    {
        $factorCallable = $this->coVarianceFactor;
        $this->factor = $factorCallable();

        // Cache covariance and inverse covariance of the data
        if (null === $this->_dataInvCov) {
            // original used the numpy.cov function.
            $this->_dataCovariance = pow(Statistics::stdev($this->dataset, true), 2);

            //$this->_dataInvCov = 1/ linalg.inv($this->_dataCovariance)
            $this->_dataInvCov = 1 / $this->_dataCovariance;
        }

        $this->covariance = $this->_dataCovariance * pow($this->factor, 2);
        $this->invCov = $this->_dataInvCov / pow($this->factor, 2);
        $this->normFactor = sqrt(2 * M_PI * $this->covariance) * count($this->dataset);
    }
929103
public static function covarianceMatrix(array $data, ?array $means = null): array
    {
        $n = count($data[0]);

        if ($means === null) {
            $means = [];
            for ($i = 0; $i < $n; ++$i) {
                $means[] = Mean::arithmetic(array_column($data, $i));
            }
        }

        $cov = [];
        for ($i = 0; $i < $n; ++$i) {
            for ($k = 0; $k < $n; ++$k) {
                if ($i > $k) {
                    $cov[$i][$k] = $cov[$k][$i];
                } else {
                    $cov[$i][$k] = self::fromDataset(
                        $data,
                        $i,
                        $k,
                        true,
                        $means[$i],
                        $means[$k]
                    );
                }
            }
        }

        return $cov;
    }
634202
def calibrateprofile():
    '''
    
    '''
    pr = profile.Profile()
    calib = []
    crepeat = 10
    for i in range(crepeat):
            calib.append(pr.calibrate(10000))
    final = sum(calib) / crepeat
    profile.Profile.bias = final # Apply computed bias to all Profile instances created hereafter
    return final
1252443
def perform_pca(A):
    
    # First subtract the mean
    M = (A-numpy.mean(A.T, axis=1)).T
    # Get eigenvectors and values of covariance matrix
    return numpy.linalg.eig(numpy.cov(M))
330547
public function covariance() : self
    {
        $b = $this->subtract($this->mean());

        return $b->matmul($b->transpose())
            ->divide($this->m);
    }
574860
def fit(self, method=None, **kwargs):
        

        cov_matrix = kwargs.get('cov_matrix', None)
        iterations = kwargs.get('iterations', 1000)
        nsims = kwargs.get('nsims', 10000)
        optimizer = kwargs.get('optimizer', 'RMSProp')
        batch_size = kwargs.get('batch_size', 12)
        mini_batch = kwargs.get('mini_batch', None)
        map_start = kwargs.get('map_start', True)
        learning_rate = kwargs.get('learning_rate', 0.001)
        record_elbo = kwargs.get('record_elbo', None)
        quiet_progress = kwargs.get('quiet_progress', False)

        if method is None:
            method = self.default_method
        elif method not in self.supported_methods:
            raise ValueError("Method not supported!")

        if method == 'MLE':
            return self._optimize_fit(self.neg_loglik, **kwargs)
        elif method == 'PML':
            return self._optimize_fit(self.neg_logposterior, **kwargs)   
        elif method == 'M-H':
            return self._mcmc_fit(nsims=nsims, method=method, cov_matrix=cov_matrix,
                map_start=map_start, quiet_progress=quiet_progress)
        elif method == "Laplace":
            return self._laplace_fit(self.neg_logposterior) 
        elif method == "BBVI":
            if mini_batch is None:
                posterior = self.neg_logposterior
            else:
                posterior = self.mb_neg_logposterior
            return self._bbvi_fit(posterior, optimizer=optimizer, iterations=iterations,
                batch_size=batch_size, mini_batch=mini_batch, map_start=map_start, 
                learning_rate=learning_rate, record_elbo=record_elbo, quiet_progress=quiet_progress)
        elif method == "OLS":
            return self._ols_fit()