Spontaneous symmetry breaking at a finite temperature $T$: How is the state dscribed as a function of $T$? Consider the equilibrium state of a statistical system with infinite DOF at a finite temperature $T$. For example, a Heisenberg ferromagnet with Hamiltonian $$H=-J\sum\limits_{i,j}\textbf{s}_i\cdot \textbf{s}_j\tag{1}$$ One can see that $\textbf{s}_i\cdot\textbf{s}_j$ is a scalar and therefore invariant under rotation. 
However, if the temperature $T$ is greater than a critical value $T_c$, then the equilibrium state respects the symmetry of the Hamiltonian. And if $T<T_c$, the symmetry of the Hamiltonian is broken by the equilibrium state. Intuitively, it is clear: a paramagnetic phase is invariant under rotations because the spins are randomly oriented and the ferromagnetic phase is not invariant under rotations because the spins align in a definite direction.
I want to make this picture mathematically precise. Can we describe the equilibrium state, in general, at temperature $T$ irrespective of whether $T>T_c$ or $T<T_c$. At finite temperature, the equilibrium state is not one of the infinite minimum energy configurations. The state is an admixture of all minimum energy states as well as all excited states weighted by appropriate Boltzmann factor.
How can I mathematically describe the equilibrium configuration so that I can explicitly see (like in Eqn. (1)) it breaks rotational invariance (symmetry of the Hamiltonian) for $T<T_c$ and restores it for $T>T_c$? 
 A: This is not a simple matter, because spontaneous symmetry breaking (and, more generally, phase transitions) can only occur for infinite systems.
I'll only discuss the classical description, and to keep things as simple as possible I'll only consider a discrete spin system on $\mathbb{Z}^d$, such as the Ising model. For such systems, the states are identified with probability measures on the set of all infinite configurations of spins.
Note that for infinite systems, you cannot write down the probability measures as being proportional to $e^{-\beta H}$, since the energy of an infinite system is in general undefined (or infinite). What always makes sense is the ratio of the probabilities of two configurations that coincide outside a finite set $\Lambda$. Indeed, the energy difference is then finite (if the interactions are absolutely summable).
This leads to the following characterization of the infinite-volume probability measures describing infinite systems of spins: $\mu$ is a Gibbs measure if it satisfies the DLR equations, that is, if 
$$
\mu(\sigma \text{ inside }\Lambda \,|\, \omega\text{ outside }\Lambda) =
\frac1{Z_{\Lambda;\beta}^\omega} e^{-\beta H(\sigma_\Lambda\omega_{\Lambda^c})}
$$
for all finite set of spins $\Lambda\subset\mathbb{Z}^d$ and (almost) all configurations $\omega$ outside $\Lambda$.
Existence of a solution of these equations is in general guaranteed. However, the interesting thing is that uniqueness fails in general. For a given Hamiltonian and at a given inverse temperature, there may be several probability measures satisfying these equations. Each of them describes one of the equilibrium states of the system.
Under very weak assumptions, one can prove that uniqueness always hold at sufficiently high temperature. But it may happen that there are several at low temperatures. This is what happens for the Ising model in dimensions $2$ and above, or for the classical Heisenberg model in dimensions $3$ and above. In the latter two cases, the symmetry of the Hamiltonian is usually broken in those Gibbs measures. For example, in the Ising model, you have two measures describing a system with positive, respectively negative, magnetization; they can be obtained by taking a finite system with $+$, respectively $-$, boundary condition and taking the thermodynamic limit. Alternatively, you can introduce a magnetic field (you have a unique state in the case) and let the intensity of the magnetic field go to $0$; if you let $h\downarrow 0$, you get the positively magnetized state, of you prefer $h\uparrow 0$, you get the negatively magnetized one.
I'll stop here, but refer you to our book on this topic, which can be downloaded here, for a much more extensive and detailed discussion.
To summarize, let me come back to your question. It is not possible to describe the state as $T$ varies (in particular as you cross the critical temperature), since spontaneous breaking of a symmetry leads to the existence of several distinct states: one for each realization of the broken symmetry, for example for each possible direction of the spontaneous magnetization. (Note that you may also get additional states.) What you can do is characterize the set of all states at a given temperature and see whether there are states in which the symmetry is broken.
Finally, similar considerations also apply to the quantum case, but with substantial additional complications.
