We present a first-principles lattice QCD calculation of the local form factors describing the $B_{s}\to \mu^{+}\mu^{-}\gamma$ decay. We focus on the region of large di-muon invariant masses $\sqrt{q^{2} }\gtrsim 4.2~{\mathrm GeV}$ , where the contributions from the four-quarks operators in the effective weak Hamiltonian (which are neglected at present) are expected to be small. We use our results for the form factors to determine the branching fraction for $B_{s} \to \mu^{+}\mu^{-}\gamma$, which has been recently measured by the LHCb in the region $\sqrt{q^{2}} > 4.9~{\rm GeV}$.
