  We extend finite-temperature tensor network methods to compute Matsubara
imaginary-time correlation functions, building on the minimally entangled
typical thermal states (METTS) and purification algorithms. While
imaginary-time correlation functions are straightforward to formulate with
these methods, care is needed to avoid convergence issues that would result
from naive estimators. As a benchmark, we study the single-band Anderson
impurity model, even though the algorithm is quite general and applies to
lattice models. The special structure of the impurity model benchmark system
and our choice of basis enable techniques such as reuse of high-probability
METTS for increasing algorithm efficiency. The results are competitive with
state-of-the-art continuous time Monte Carlo. We discuss the behavior of
computation time and error as a function of the number of purified sites in the
Hamiltonian.
