Population studies of -ray blazars were undertaken soon after the recognition of the -ray blazar class with EGRET .  performed a < V / Vmax > analysis assuming no density evolution and showed that luminosity evolution of EGRET blazars was implied by the data. With a larger data set, and using radio data to ensure the sample was unbiased in regard to redshift determination,  again found that luminosity evolution was required. They obtained best-fit values through the maximum likelihood method that gave an AGN contribution to the EGRET -ray background at the level of 25%.
 postulated a radio / -ray correlation in blazars, and tried to correct for the duty cycle and -ray spectral hardening of flaring states. They found that essentially 100% of the EGRET -ray background arises from unresolved blazars and AGNs. In later work , they predict that GLAST will detect 5000 blazars to a flux level of 2 × 10-9 ph(> 100 MeV)/(cm2-s), which will be reached with GLAST after 4 years. They did not, however, fit the blazar redshift distribution to provide a check on their model, nor distinguish between flat spectrum radio quasar (FSRQ) and BL Lac objects.
The crucial underlying assumption of this approach, which has been developed in recent work [18, 33], is that there is a simple relation between the radio and -ray fluxes of blazars. Because a large number of EGRET -ray blazars (primarily FSRQs) are found in the 5 GHz, > 1 Jy  catalog, a radio / -ray correlation is expected. This correlation is not, however, evident in 2.7 and 5 GHz monitoring of EGRET -ray blazars . X-ray selected BL objects are also not well-sampled in GHz radio surveys. Studies based on correlations between the radio and -ray emissions from blazars must therefore consider the very different properties and histories of FSRQs and BLs and their separate contributions to the -ray background.
Treatments of blazar statistics that avoid any radio / -ray correlation and separately consider FSRQs and BL Lac objects have been developed by  and . In the  study, blazar spectra were calculated assuming an injection electron number index of -2. Distributions in injected particle energy in BL Lac and FSRQ jets were separately considered, with a simple description of density evolution given in the form of a cutoff at some maximum redshift zmax. Depending on the value of zmax,  concluded that as much as 40 - 80% of the EGRB is produced by unresolved AGNs, with 70 - 90% of the emission from FR 1 galaxies and BL Lac objects.
In my recent study , I also use a physical model to fit the EGRET data on the redshift and size distribution of EGRET blazars. The EGRET blazar sample consists of 46 FSRQs and 14 BL Lac objects that were detected in the Phase 1 EGRET all-sky survey , with fluxes as reported in the Third EGRET catalog . A blazar is approximated by a relativistic spherical ball entraining a tangled magnetic field and containing an isotropic, power-law distribution of nonthermal electrons. Single electron power-law distributions were used in the study, with indices p = 3.4 for FSRQs and p = 3.0 for BL Lac objects, giving spectral indices = -0.2 and = 0.0, respectively, as shown by observations [31, 50]. Beaming patterns appropriate to external Compton and synchrotron self-Compton processes, and bulk Lorentz factor = 10 and = 4, were used in FSRQs and BL Lac objects, respectively. The comoving directional luminosities l'e and blazar comoving rate densities (blazar formation rate; BFRs) for the two classes were adjusted to give agreement with the data. The threshold detector sensitivity -8, in units of 10-8 ph(> 100 MeV)/(cm2-s), was nominally taken to be -8 = 15 for the two-week on-axis EGRET sensitivity, and -8 = 0.4 for the one-year all-sky sensitivity of GLAST. Due to incompleteness of the sample near threshold, the EGRET threshold was adjusted to -8 = 25. Because a mono-luminosity function was used, the range in apparent powers is entirely kinematic in this model, arising from the different, randomly oriented jet directions.
By using a minimalist blazar model, the model parameters were severely constrained. The FSRQ data were fit with l'e = 1040 ergs/(s-sr) and a BFR that was 15 &215; greater at z 2 - 3 than at present. The BL Lac data, by contrast, could not be fit using a fixed luminosity. A model that could jointly fit the redshift and size distribution of BL Lac objects required that BL Lac objects be brighter and less numerous that in the past, consistent with a picture where FSRQs evolve into BL Lac objects [8, 27].
Fig. 1b shows the fitted EGRET redshift distributions and predicted redshift distributions of galaxies and blazars at different GLAST sensitivities . The fits to the EGRET size distributions of FSRQs and BL Lac objects, and extrapolations of the model size distributions to lower flux thresholds, are shown in Fig. 1c. After one year of observations with GLAST (-8 0.4), 800 FSRQs/FR2 and 200 BL Lac/FR1 galaxies and -ray blazars are predicted. This is a lower prediction, and additional hard-spectrum blazars to which EGRET was not sensitive could increase this number, but not by more than a factor 2. The contribution of unresolved blazars below a flux level of -8 12.5 - 25 to the EGRB is shown in Fig. 1a. As can be seen, the total blazar / galaxy contribution is less than 20 - 30% of the EGRET EGRB intensity, meaning that other classes of sources must make a significant contribution.