In this article, a mathematical model has been developed to predict the non-linear dynamic behaviour of a high-speed unbalanced rotating shaft due to ball size variations and varying number of balls. In the mathematical formulation, the contact between balls and inner/outer races is considered as a non-linear spring, whose stiffness is obtained using Hertzian elastic deformation theory. A detailed contact-damping model reflecting the influences of the surface profiles and the speeds of both contacting elements is developed and applied in the rolling element bearing formulation. This formulation also accounts for sources of non-linearity such as Hertzian contact force, varying number of balls, and off-sized balls resulting transition from noncontact to contact state between balls and races. The equations of motion of a rolling element bearing are formulated using Lagrange's equation considering the vibration characteristics of individual constituents such as inner race, outer race, rolling elements, and shaft. To overcome the strong non-linearity of the governing equations of motion, a modified Newmark-β time integration technique has been used to solve the equations of motion numerically. All results have been presented in the form of fast Fourier transformations and Poincaré maps. The highest radial vibrations due to ball size variation are at a speed of the number of balls times the cage speed (ω = kω cage Hz). The other vibrations due to ball size variation also occurr at (X ± kω cage ), where k is a constant. When the number of balls is increased, the centre of oscillations approaches zero, implying a stiffer system. From this it can be predicted that increasing the number of balls will reduce the effect of the BPF. Hence it is implied from the analysis that the number of balls and the size variation of balls are important parameters and should be considered at the design stage.