How can holonomy corrections be introduced in $f(R)$ gravity?

We study the introduction of holonomy corrections in $f(R)$ gravity. We will show that there are infinitely many ways, as many as canonical transformations, to introduce this kind of corrections, depending on the canonical variables (two coordinates and its conjugate momenta) used to obtain the Hamiltonian. In each case, these corrections lead, at effective level, to different modified holonomy corrected Friedmann equations in $f(R)$ gravity, which are in practice analytically unworkable, i.e. only numerical analysis can be used to understand its dynamics. Finally, we give arguments in favour of one preferred set of variables, the one that conformally maps $f(R)$ to Einstein gravity, because for these variables the dynamics of the system has a clear physical meaning: the same as in standard Loop Quantum Cosmology, where the effective dynamics of a system can be analytically studied.