We develop in this paper efficient and robust numerical methods for solving anisotropic Cahn-Hilliard systems. We construct energy stable schemes for the time discretization of the highly nonlinear anisotropic Cahn-Hilliard systems by using a stabilization technique. At each time step, these schemes lead to a sequence of linear coupled elliptic equations with constant coefficients that can be efficiently solved by using a spectral-Galerkin method. We present numerical results that are consistent with earlier work on this topic, and also carry out various simulations, such as the linear bi-Laplacian regularization and the nonlinear Willmore regularization, to demonstrate the efficiency and robustness of the new schemes.