In the study of stochastic PDEs with colored, non-trace class space-time noise, one frequently encounters Gaussian series of the form $$g \sum_{n\geq 1} γ_n μ_n f_n, $$ where $(γ_n)_{n}$ is a sequence of standard independent Gaussian variables, $g$ is an $L^η(\mathcal{O})$ function, $(μ_n)_{n}$ is a sequence of scalars, and $(f_n)_n$ is an orthonormal system in $L^2(\mathcal{O})$ where $\mathcal{O} \subseteq \mathbb{R}^d$ is an open set. In this manuscript, we establish necessary and sufficient conditions for the above sum to converge in Bessel potential spaces $H^{-s,q}(\mathcal{O})$. The latter can be interpreted as a Sobolev embedding for Gaussian series. Our main theorem is formulated using weighted sequence spaces that encode the $L^\infty$-growth of the orthonormal system $(f_n)_{n}$, a feature that is crucial for obtaining sharp estimates. We apply our results to the stochastic heat equation with additive non-trace class noise. In this case, our conditions capture the scaling relationship between the heat operator and the coloring of the noise.