We study universality properties of the Epstein zeta function En(L, s) for lattices L of large dimension n and suitable regions of complex numbers s. Our main result is that, as n → ∞, En(L, s) is universal in the right half of the critical strip as L varies over all n-dimensional lattices L. The proof uses an approximation result for Dirichlet polynomials together with a recent result on the distribution of lengths of lattice vectors in a random lattice of large dimension and a strong uniform estimate for the error term in the generalized circle problem. Using the same approach we also prove that, as n → ∞, En(L1, s) − En(L2, s) is universal in the full half-plane to the right of the critical line as (L1, L2) varies over all pairs of n-dimensional lattices. Finally, we prove a more classical universality result for En(L, s) in the s-variable valid for almost all lattices L of dimension n. As part of the proof we obtain a strong bound of En(L, s) on the critical line that is subconvex for n ≥ 5 and almost all n-dimensional lattices L.Let us also note that, since ζ(s) has an Euler product, it is necessary to assume that the function f in Theorem 1.1 is nonvanishing on K. However, for zeta functions