In this work, based on first order shear deformation theory, the equations of motion for a viscoelastic rectangular plate have been derived, which consider the effects of shear deformation and rotary inertia. For the viscoelastic material, the Kelvin-Voigt model was used. Then, using the Galerkin method, expanding the invariant manifold theory for plate problems, and using the multiple scale method, the closed form relations for amplitude, non-linear, natural frequencies, internal resonances, and non-linear mode shapes have been obtained. Finally, as an example, a plate with simply supported boundary condition has been considered and using the obtained relations, the non-linear natural frequency and mode shape have been calculated. Then, the effects of different system parameters, such as thickness, Poisson's ratio, and damping coefficient on the response of the system have been investigated.