In this study, an efficient polygonal finite element method (PFEM) in combination with quadratic serendipity shape functions is proposed to study nonlinear static and dynamic responses of functionally graded (FG) plates with porosities. Two different porosity types including even and uneven distributions through the plate thickness are considered. The quadratic serendipity shape functions over arbitrary polygonal elements including triangular and quadrilateral ones, which are constructed based on a pairwise product of linear shape functions, are employed to interpolate the bending strains. Meanwhile, the shear strains are defined according to the Wachspress coordinates. By using the Timoshenko's beam to interpolate the assumption of the strain field along the edges of polygonal element, the shear locking phenomenon can be naturally eliminated. Furthermore, the C0–type higher-order shear deformation theory (C0–HSDT), in which two additional variables are included in the displacement field, significantly improves the accuracy of numerical results. The nonlinear equations of static and dynamic problems are solved by Newton–Raphson iterative procedure and by Newmark's integration scheme in association with the Picard methods, respectively. Through various numerical examples in which complex geometries and different boundary conditions are involved, the proposed approach yields more stable and accurate results than those generated using other existing approaches.