Presented is the search for anomalous Higgs boson decays into two axion-like particles (ALPs) using the full Run 2 data set of 140 fb$^{-1}$ of proton-proton collisions at a centre-of-mass energy of 13 TeV recorded by the ATLAS experiment. The ALPs are assumed to decay into two photons, providing sensitivity to recently proposed models that could explain the $(g-2)_\mu$ discrepancy. This analysis covers an ALP mass range from 100 MeV to 62 GeV and ALP-photon couplings in the range $10^{-5}\, \text{TeV}^{-1}<C_{a\gamma\gamma}/\Lambda<1\, \text{TeV}^{-1}$, and therefore includes signatures with significantly displaced vertices and highly collinear photons. No significant excess of events above the Standard Model background is observed. Upper limits at 95% confidence level are placed on the branching ratio of the Higgs boson to two ALPs in the four-photon final state, and are in the range of $ 10^{-5}$ to $3\times 10^{-2}$, depending on the hypothesized ALP mass and ALP-photon coupling strength.