Calculate the Constant-Correlation covariance matrix.

        :return: shrunk sample covariance matrix
        :rtype: np.ndarray

def shrink(self):
        """
        Calculate the Constant-Correlation covariance matrix.

        :return: shrunk sample covariance matrix
        :rtype: np.ndarray
        """
        x = np.nan_to_num(self.X.values)

        # de-mean returns
        t, n = np.shape(x)
        meanx = x.mean(axis=0)
        x = x - np.tile(meanx, (t, 1))

        # compute sample covariance matrix
        sample = (1.0 / t) * np.dot(x.T, x)

        # compute prior
        var = np.diag(sample).reshape(-1, 1)
        sqrtvar = np.sqrt(var)
        _var = np.tile(var, (n,))
        _sqrtvar = np.tile(sqrtvar, (n,))
        r_bar = (np.sum(sample / (_sqrtvar * _sqrtvar.T)) - n) / (n * (n - 1))
        prior = r_bar * (_sqrtvar * _sqrtvar.T)
        prior[np.eye(n) == 1] = var.reshape(-1)

        # compute shrinkage parameters and constant
        if self.delta is None:

            # what we call pi-hat
            y = x ** 2.0
            phi_mat = np.dot(y.T, y) / t - 2 * np.dot(x.T, x) * sample / t + sample ** 2
            phi = np.sum(phi_mat)

            # what we call rho-hat
            term1 = np.dot((x ** 3).T, x) / t
            help_ = np.dot(x.T, x) / t
            help_diag = np.diag(help_)
            term2 = np.tile(help_diag, (n, 1)).T * sample
            term3 = help_ * _var
            term4 = _var * sample
            theta_mat = term1 - term2 - term3 + term4
            theta_mat[np.eye(n) == 1] = np.zeros(n)
            rho = sum(np.diag(phi_mat)) + r_bar * np.sum(
                np.dot((1.0 / sqrtvar), sqrtvar.T) * theta_mat
            )

            # what we call gamma-hat
            gamma = np.linalg.norm(sample - prior, "fro") ** 2

            # compute shrinkage constant
            kappa = (phi - rho) / gamma
            shrinkage = max(0.0, min(1.0, kappa / t))
            self.delta = shrinkage
        else:
            # use specified constant
            shrinkage = self.delta

        # compute the estimator
        sigma = shrinkage * prior + (1 - shrinkage) * sample
        return self.format_and_annualise(sigma)