We consider centered conditionally Gaussian d-dimensional vectors X with random covariance matrix Ξ having an arbitrary probability distribution law on the set of nonnegative definite symmetric d × d matrices Md+. The paper deals with the evaluation problem of mean values \( E\left[ {\prod\nolimits_{i = 1}^{2n} {\left( {{c_i},X} \right)} } \right] \) for ci ∈ ℝd, i = 1, …, 2n, extending the Wick theorem for a wide class of non-Gaussian distributions. We discuss in more detail the cases where the probability law ℒ(Ξ) is infinitely divisible, the Wishart distribution, or the inverse Wishart distribution. An example with Ξ\( = \sum\nolimits_{j = 1}^m {{Z_j}{\sum_j}} \), where random variables Zj, j = 1, …, m, are nonnegative, and Σj ∈ Md+, j = 1, …, m, are fixed, includes recent results from Vignat and Bhatnagar, 2008.