Predicting coupled frequency response of piezo-elasto-magnetic (PEM) structures is crucial for accurately designing sensors/actuators, energy harvesters and other smart structures. The numerous coupling parameters involved make the numerical analysis through the conventional finite element methods (FEM) cumbersome and time-consuming, particularly for non-uniformly shaped/variable thickness structures. Hence, in this article, a hybrid approach integrating the computational benefits of FEM and artificial neural network (ANN) models has been proposed to predict the coupled nonlinear frequency response (NLFR) of porous functionally graded PEM (PFG-PEM) plates with non-uniform geometries. Through this approach, the computational efforts are substantially reduced retaining appreciable accuracy. A FEM model based on Hamilton's principle, von Karman's nonlinearity and higher-order shear deformation theory (HSDT) was initially developed for non-uniform PFG-PEM plates. The large datasets collected from the nonlinear FEM simulation are used to train an ANN model that can accurately predict the NLFR of non-uniform PFG-PEM plates for out-of-range input data sets. The plates in the current study have non-uniform thicknesses varying bi-linearly, linearly, and exponentially. The different variants of PFG-PEM composites and porosity patterns are evaluated, whose material property varies across the thickness according to a power law distribution. In addition, two forms of electromagnetic boundary conditions, such as open and closed circuits, are enforced on the plate, and its NLFR is assessed. Further, several numerical examples are presented to understand the interdependency of several material and geometrical parameters on the overall NLFR of PFG-PEM plates. This predictive tool can be readily used for further optimisation of smart structural design, significantly reducing the time consumed.