  We present a general theory of classical metastability in open quantum
systems. Metastability is a consequence of a large separation in timescales in
the dynamics, leading to the existence of a regime when states of the system
appear stationary, before eventual relaxation toward a true stationary state at
much larger times. In this work, we focus on the emergence of classical
metastability, i.e., when metastable states of an open quantum system with
separation of timescales can be approximated as probabilistic mixtures of a
finite number of states. We find that a number of classical features follow
from this approximation, for the manifold of metastable states, long-time
dynamics between them, and symmetries of the dynamics. Namely, those states are
approximately disjoint and thus play the role of metastable phases, the
relaxation toward the stationary state is approximated by a classical
stochastic dynamics between them, and weak symmetries correspond to their
permutations. Importantly, the classical dynamics is observed not only on
average, but also at the level of individual quantum trajectories: We show that
time coarse-grained continuous measurement records can be viewed as noisy
classical trajectories, while their statistics can be approximated by that of
the classical dynamics. Among others, this explains how first-order dynamical
phase transitions arise from metastability. Finally, to verify the presence of
classical metastability in a given open quantum system, we develop an efficient
numerical approach that delivers the set of metastable phases together with the
effective classical dynamics. Since the proximity to a first-order dissipative
phase transition manifests as metastability, the theory and tools introduced in
this work can be used to investigate such transitions through the metastable
behavior of many-body systems of moderate sizes accessible to numerics.
