  We demonstrate that a wide range of viable f(R) parameterizations (including
the Hu & Sawicki and the Starobinsky models) can be expressed as perturbations
deviating from the LCDM Lagrangian. We constrain the deviation parameter b
using a combination of geometrical and dynamical observational probes. In
particular, we perform a joint likelihood analysis of the recent Supernovae
Type Ia data, the Cosmic Microwave Background shift parameters, the Baryonic
Acoustic Oscillations and the growth rate data provided by the various galaxy
surveys. This analysis provides constraints for the following parameters: the
matter density Omega_{m0}, the deviation from LCDM parameter b and the growth
index gamma(z). We parametrize the growth index gamma(z) in three manners
(constant, Taylor expansion around z=0, and Taylor expansion around the scale
factor). We point out the numerical difficulty for solving the generalized f(R)
Friedman equation at high redshifts due to stiffness of the resulting ordinary
differential equation. We resolve this problem by constructing an efficient
analytical perturbative method in the deviation parameter b. We demonstrate
that this method is highly accurate, by comparing the resulting analytical
expressions for the Hubble parameter, with the numerical solutions at low and
intermediate redshifts. Surprisingly, despite of its perturbative nature, the
accuracy of the method persists even for values of b that are of O(1).
