Glueballs, predicted by QCD, are states predominantly made of gluons, yet no clear experimental observation has been made since their prediction. Recently, the BESIII collaboration reported the discovery of a pseudoscalar glueball-like particle X(2370), compatible with lattice QCD predictions. Lattice QCD is ideal for investigating these composite states, but challenges such as the signal-to-noise ratio and the analysis with multi-meson interpolating operators hinder comprehensive studies of the glueball spectrum in full QCD. I will discuss our recent use of the multilevel algorithm to improve the signal-to-noise ratio in pure Yang-Mills theory simulations. The multilevel algorithm enables us to obtain signal at larger distances (up to $\sim$ 0.8 fm) compared to standard algorithms at the same cost, which will be crucial for studies of glueball-charmonium mixing in full QCD.