We investigate particle spectra and elliptic flow coefficients in relativistic heavy ion collisions by taking into account the distortion of phase space distributions due to bulk viscosity at freezeout. We first calculate the distortion of phase space distributions in a multi-component system within the Grad's fourteen moment method. We find some subtle issues when one matches macroscopic variables with microscopic momentum distributions in a multi-component system and develop a consistent procedure to uniquely determine the corrections in the phase space distributions. Next, we calculate particle spectra by using the Cooper-Frye formula to see the effect of the bulk viscosity. In spite of the relative smallness of the bulk viscosity, we find that it is likely to have a visible effect in particle spectra and elliptic flow coefficients. This indicates the importance of bulk viscosity together with shear viscosity if one wants to constrain the transport coefficients with better accuracy from experimental data.