Understanding how the properties of heavy mesons change as temperature increases is crucial for gaining valuable insights into the quark-gluon plasma. Information about meson masses and decay widths is encoded in the meson spectral function, which, in principle, can be extracted from Euclidean correlation functions via generalised Laplace transformations. However, this inverse problem is ill-posed for lattice correlation functions and requires regularisation. In this work, we present the latest results for bottomonium spectral functions obtained within the lattice NRQCD framework using the Backus-Gilbert regularisation, along with two other variants, one of which is commonly referred to as the HLT method. Our analysis employs Generation 2L anisotropic lattice configurations produced by the \textsc{Fastsum} collaboration.