The effects of atmospheric icing can be anticipated by Computational Fluid Dynamics (CFD). Past studies show that the convective heat transfer influences the ice accretion and is itself a function of surface roughness. Uncertainty quantification (UQ) could help quantify the impact of surface roughness parameters on the reliability of ice accretion prediction. This paper aims to quantify ice accretion uncertainties and identify the key surface roughness correction parameters contributing the most to the uncertainties in a Reynolds-Averaged Navier-Stokes (RANS) formulation. Ice accretion simulations over a rough flat plate using two thermal correction models are used to construct a RANS database. Non-Intrusive Polynomial Chaos Expansion (NIPCE) metamodels are developed to predict the convective heat transfer and icing characteristics of the RANS database. The metamodels allow for the computation of the 95% confidence intervals of the output probability distribution (PDF) and of the sensitivity indexes of the roughness parameters according to their level of influence on the outputs. For one of the thermal correction models, the most influential parameter is the roughness height, whereas for the second model it is the surface correction coefficient. In addition, the uncertainty on the freestream temperature has a minor impact on the ice accretion sensitivity compared to the uncertainty on the roughness parameters.