We apply non-perturbative renormalization to bilinears composed of improved staggered fermions. We explain how to generalize the method to staggered fermions in a way which is consistent with the lattice symmetries, and introduce a new type of lattic
e bilinear which transforms covariantly and avoids mixing. We derive the consequences of lattice symmetries for the propagator and vertices. We implement the method numerically for hypercubic-smeared (HYP) and asqtad valence fermion actions, using lattices with asqtad sea quarks generated by the MILC collaboration. We compare the non-perturbative results so obtained to those from perturbation theory, using both scale-independent ratios of bilinears (of which we calculate 26), and the scale-dependent bilinears themselves. Overall, we find that one-loop perturbation theory provides a successful description of the results for HYP-fermions if we allow for a truncation error of roughly the size of the square of the one-loop term (for ratios) or of size O(1) times alpha^2 (for the bilinears themselves). Perturbation theory is, however, less successful at describing the non-perturbative asqtad results.
