Is there an analogy for Wilson loops/lines in statistical mechanics?

Gauge theories often appear in condensed matter. For introduction I recommend read brilliant introduction article by professor Wen Topological order: from long-range entangled quantum matter to an unification of light and electrons.

In gauge theory, observables are gauge-invariant operators (because gauge symmetry is not symmetry! It is redundancy in our description.). In gauge theory without matter such objects are Wilson and t'Hooft loops. So in any theory with gauge symmetry this observables are natural and very important.

For concrete examples consult with Field Theories of Condensed Matter Physics. See 8.4 about disordered spin states, 9.6 the Ising gauge theory, 9.9 Boundary conditions and topology,....

Ising gauge theory

Here I follow Z 2 gauge theory by Subir Sachdev.

System is defined through Hamiltonian:

There are two gapped phases in the theory, which are necessarily separated by a phase transition. Remarkably, unlike all previously known cases, this phase transition was not required by the presence of a broken symmetry in one of the phases: there was no local order parameter characterizing the phase transition. Instead, Wegner argued for the presence of a phase transition using the behavior of the Wegner-Wilson loop operator $W_C$, which is the product of $\sigma_z$ on the links of any closed contour C on the direct square lattice:

In this example, Wilson loop is nonlocal order parameter for phase transition.

Suppose you have a system with physical field $\sigma(\vec{r})$, which is constrained by some law $L[\sigma](\vec{r})=0$. Naively, to do statistical mechanics, one has to carry around the constraint $L[\sigma]=0$ in all computations.

In some cases, in particular when $L$ is linear and built out of derivatives, the constraint can be explicitly solved by writing $\sigma = K[\phi]$, in terms of a new field $\phi$. This allows us to work directly on the constrained manifold. The price of this representation is that this new field is, in general, a gauge field. This gives a general mechanism for gauge theory in statistical mechanics.

For example, if $\sigma=\sigma_{ij}$ is the stress tensor then the condition of mechanical equilibrium is $\partial_i \sigma_{ij}=0$. This can be solved (in 3 dimensions) by writing $\sigma_{ij} = \epsilon_{ikl} \partial_k \phi_{lj}$. The new field $\phi$ cannot directly be a physical field, since we can add any function $\Delta \phi_{lj}$ for which $\epsilon_{ikl} \partial_k \Delta \phi_{lj}=0$, such as constants. We are thus led to consider the gauge theory for $\phi$. The original physical Lagrangian, once written in terms of $\phi$, will have Wilson-loop-like terms, after using Stokes' theorem. For example, the surface integral $\int_{A} dS n_i \sigma_{ij} = \int_{\partial A} dl_i \phi_{ij}$, where $A$ is some surface with boundary $\partial A$.

The example here is taken from the physics of amorphous solids at low temperature, for which the states of mechanical equilibrium control the vibrational and rheological properties of the solid. In this case one needs to also impose torque balance, so $\phi$ is further rewritten in terms of a more primitive field, but the result is the same: gauge theory.