We compare the color-evaporation model (CEM) and nonrelativistic QCD (NRQCD) factorization predictions for inclusive quarkonium production. Using the NRQCD factorization formulas for quarkonium production and for perturbative QQ production, we deduce relationships that are implied by the CEM between the nonperturbative NRQCD matrix elements that appear in the factorization formula for quarkonium production. These relationships are at odds with the phenomenological values of the matrix elements that have been extracted from the Tevatron data for charmonium production at large transverse momentum. A direct comparison of the CEM and NRQCD factorization predictions with the CDF charmonium production data reveals that the CEM fits to the data are generally unsatisfactory, while the NRQCD factorization fits are generally compatible with the data. The inclusion of k T smearing improves the CEM fits substantially, but significant incompatibilities remain. The NRQCD factorization fits to the χ c data indicate that multiple gluon radiation is an essential ingredient in obtaining the correct shape of the cross section as a function of p T .