  We study the distortions of equilibrium spectra of relic neutrinos due to the
interactions with electrons, positrons, and neutrinos in the early Universe. We
solve the integro-differential kinetic equations for the neutrino density
matrix, including three-flavor oscillations and finite temperature corrections
from QED up to the next-to-leading order $\mathcal{O}(e^3)$ for the first time.
In addition, the equivalent kinetic equations in the mass basis of neutrinos
are directly solved, and we numerically evaluate the distortions of the
neutrino spectra in the mass basis as well, which can be easily extrapolated
into those for non-relativistic neutrinos in the current Universe. In both
bases, we find the same value of the effective number of neutrinos, $N_{\rm
eff} = 3.044$, which parameterizes the total neutrino energy density. The
estimated error for the value of $N_{\rm eff}$ due to the numerical
calculations and the choice of neutrino mixing parameters would be at most
0.0005.
