  We address the cosmological role of an additional ${\cal O}(1)$ eV sterile
neutrino in modified gravity models. We confront the present cosmological data
with predictions of the FLRW cosmological model based on a variant of $f(R)$
modified gravity proposed by one of the authors previously. This viable
cosmological model which deviation from general relativity with a cosmological
constant $\Lambda$ decreases as $R^{-2n}$ for large, but not too large values
of the Ricci scalar $R$ provides an alternative explanation of present dark
energy and the accelerated expansion of the Universe. Various up-to-date
cosmological data sets exploited include Planck CMB anisotropy, CMB lensing
potential, BAO, cluster mass function and Hubble constant measurements. We find
that the CMB+BAO constraints strongly the sum of neutrino masses from above.
This excludes values $\lambda\sim 1$ for which distinctive cosmological
features of the model are mostly pronounced as compared to the $\Lambda$CDM
model, since then free streaming damping of perturbations due to neutrino rest
masses is not sufficient to compensate their extra growth occurring in $f(R)$
gravity. Thus, we obtain $\lambda>8.2$ ($2\sigma$) with cluster systematics and
$\lambda>9.4$ ($2\sigma$) without that. In the latter case we find for the
sterile neutrino mass
$0.47\,\,\rm{eV}$$\,<\,$$m_{\nu,\,\rm{sterile}}$$\,<\,$$1\,\,\rm{eV}$
($2\sigma$) assuming the active neutrinos are massless, not significantly
larger than in the standard $\Lambda$CDM with the same data set:
$0.45\,\,\rm{eV}$$\,<\,$$m_{\nu,\,\rm{sterile}}$$\,<\,$$0.92\,\,\rm{eV}$
($2\sigma$). However, a possible discovery of a sterile neutrino with the mass
$m_{\nu,\,\rm{sterile}} \approx 1.5\,$eV motivated by various anomalies in
neutrino oscillation experiments would favor cosmology based on $f(R)$ gravity
rather than the $\Lambda$CDM model.
