In this paper, based on the L2-1Ï scheme and nonconforming EQ1rot finite element method (FEM), a numerical approximation is presented for a class of two-dimensional time-fractional diffusion equations involving variable coefficients. A novel and detailed analysis of the equations with an initial singularity is described on anisotropic meshes. The fully discrete scheme is shown to be unconditionally stable, and optimal second-order accuracy for convergence and superconvergence can be achieved in both time and space directions. Finally, the obtained numerical results are compared with the theoretical analysis, which verifies the accuracy of the proposed method.