'Estimate the full Wishart distribution parameters.
Parameters
X : array-like, shape (n_samples, n_features)
nk : array-like, shape (n_components,)
xk : array-like, shape (n_components, n_features)
sk : array-like, shape (n_components, n_features, n_features)'

def _estimate_wishart_full(self, nk, xk, sk):

    (_, n_features) = xk.shape
    self.degrees_of_freedom_ = (self.degrees_of_freedom_prior_ + nk)
    self.covariances_ = np.empty((self.n_components, n_features, n_features))
    for k in range(self.n_components):
        diff = (xk[k] - self.mean_prior_)
        self.covariances_[k] = ((self.covariance_prior_ + (nk[k] * sk[k])) + (((nk[k] * self.mean_precision_prior_) / self.mean_precision_[k]) * np.outer(diff, diff)))
    self.covariances_ /= self.degrees_of_freedom_[:, np.newaxis, np.newaxis]
