  We study black hole linear perturbation theory in a four-dimensional
Schwarzschild (anti) de Sitter background. When dealing with a positive
cosmological constant, the corresponding spectral problem is solved
systematically via the Nekrasov-Shatashvili functions or, equivalently,
classical Virasoro conformal blocks. However, this approach can be more
complicated to implement for certain perturbations if the cosmological constant
is negative. For these cases, we propose an alternative method to set up
perturbation theory for both small and large black holes in an analytical
manner. Our analysis reveals a new underlying recursive structure that involves
multiple polylogarithms. We focus on gravitational, electromagnetic, and
conformally coupled scalar perturbations subject to Dirichlet and Robin
boundary conditions. The low-lying modes of the scalar sector of gravitational
perturbations and its hydrodynamic limit are studied in detail.
