Background
No SSB in finite systems
Consider a system interacting with a heat bath at inverse temperature $\beta$, with the resultant dynamics of the system described by a Liouvillian superoperator $\mathcal{L}$. If this system is finite, then under fairly general conditions on $\mathcal{L}$, we expect the equilibrium state $\rho$, meaning $\mathcal{L}(\rho) = 0$, to be uniquely given by the Gibbs state
$$ \rho_{\mathrm{gibbs}} = \dfrac{e^{-\beta H}}{\mathrm{Tr}\left[ e^{-\beta H} \right]}, $$
where $H$ is the Hamiltonian of the system.
Let $\mathcal{G}$ be the symmetry group of the Hamiltonian, meaning there is a unitary representation $U$ of $\mathcal{G}$ such that $[H, U(g)] = 0$ for all $g \in \mathcal{G}$. It is clear that we also have $[\rho_{\mathrm{gibbs}}, U(g)] = 0$, so the Gibbs state preserves all the symmetries of the Hamiltonian. In this sense, it seems that there cannot be spontaneous symmetry breaking (SSB) in finite systems.
SSB in infinite systems (and KMS states)
The typical narrative then proceeds to say that, in fact, SSB can occur, but only in infinite systems. Here there is no guarantee that $e^{-\beta H}$ is trace-class, so in general the Gibbs state is not well defined. To extend the notion of a "thermal" state to infinite systems, one usually defines the so-called KMS states. These are the states $\phi$ which satisfy the KMS condition, which can (informally) be stated as
$$ \langle A (t) B \rangle_{\phi} = \langle B(t + i\beta) A \rangle_{\phi}, $$
for all operators $A$ and $B$ in the operator algebra, where $\langle \cdot \rangle_{\phi}$ indicates an expectation value with respect to the state $\phi$. (I omit all $C^{*}$-algebraic details here for brevity.)
There is a large body of literature showing that KMS states preserve the properties that we consider key to the definition of a thermal state, such as being equilibrium states, but remain well-defined for infinite systems.
For finite systems, I believe the KMS condition uniquely specifies a state: the Gibbs state. However, for infinite systems this is not necessarily the case, and, roughly speaking, SSB occurs when there are multiple KMS states, each of which is not preserved by the symmetry group of the Hamiltonian.
Question
Both experiments and numerical simulations show systems with behaviour that seems very similar to that of SSB (ferromagnets exist!). However, these real-world systems are clearly finite, so the above arguments would suggest that they cannot truly display SSB. What is the explanation for this discrepancy?
Thoughts on an answer
Though finite, real-world experiments can often by fairly effectively described by taking the infinite size limit. If this is appropriate, then perhaps the dynamics of these large finite systems can be well approximated by infinite systems, at least up to some large timescale $\tau$ which presumably grows quickly with system size. Then we might expect these finite systems to display signatures of SSB over the timescale $\tau$, after which they will decay to the Gibbs state and the symmetry will be restored. If this is along the right lines, can any of this be made precise?
