# Sample from parameters' posteriori distribution

Assume a linear model for the data : $$y = x^T\beta + \epsilon$$, with $$\epsilon \sim \mathcal{N}(0,\sigma^2)$$ and an a priori on the parameter: $$\beta \sim \mathcal{N}(0,\Sigma^2 _p)$$. Let $$\mathbf{y}, \mathbf{X}$$ be the training set. It is widely known that the a posteriori distribution of $$\beta$$ is $$p(\beta |\mathbf{y},\mathbf{X}) = \mathcal{N}(\sigma^{-2}(\sigma^{-2}XX^T + \Sigma_p)^{-1}Xy, \sigma^{-2} XX^T + \Sigma_p^{-1})$$

now what if we didn't know how to derive this result, how can we 'observe' this a posteriori by sampling?