Solving linear matrix equations has various applications in control theory, in engineering, in scientific computations and various other fields. By applying generalization of the Hermitian and skew-Hermitian splitting (GHSS) iteration and the hierarchical identification principle, we propose a gradient-based iterative method for finding the solution of the generalized coupled Sylvester matrix equations (including (coupled) Sylvester and Lyapunov matrix equations as special cases). We prove that the iterative solution consistently converges to the solution for any initial matrix. Some numerical examples and applications are provided to illustrate the effectiveness of the method.