Continuous change in wrinkle patterns process of thin films to a gradient substrate is the most challenging problem regarding applying reliable and robust numerical methods for postbuckling analysis of the film/substrate system. For example, in the finite element method, the postbuckling simulation suffers from the convergence issue, while, in spectral methods, it is very difficult to capture the localized behavior in soft matters when the boundary conditions are complex. When a thin film is compressed, it can form a wrinkle of a certain amount of wavelength when the compression exceeds a critical value. The compressed compliant substrate system translates to sinusoidal wrinkles and then to period-doubling wrinkling after further compression. In this work, we investigate the mathematical model arising from the changing nature of wrinkle patterns of postbuckled thin films using a robust and efficient numerical algorithm based on the spectral method to evolve the wrinkle patterns. We consider the gradient substrates of three typical variations in the modulus, namely, the symmetry, exponential, and power-law model. It has been observed that the stable equilibrium path has two bifurcation points. At the first bifurcation point, the buckling instability of wrinkling occurs, while the period-doubling buckling instability occurs at the second bifurcation point. For the substrates of material gradients of various types, the amplitude and wavelength are obtained. This study may help in better understanding of wrinkle patterns formation which could be very useful for the designing of stretchable and flexible electronic devices of most substrate systems and to avoid resonance in the noise environment.