This study investigates the numerical stability of the two-dimensional lattice Boltzmann (LB) scheme for nonlinear flows. While the BGK (single relaxation time) scheme is known for its simplicity, it can be unstable for complex flows. Conversely, the LB scheme with multiple relaxation times (MRT) offers greater stability for such nonlinear problems. However, traditional stability analysis methods like von Neumann local analysis are limited to the linear case. Here, we examine the stability of a specific nonlinear test case with fixed viscosity. We employ a decision tree, a machine learning technique valued for its interpretability, to explore and characterize the stability zone for the free relaxation parameters. To further investigate, a simple global optimization method (evolutionary algorithm) is used to identify a set of stable relaxation parameters for various test cases, including the doubly periodic shear layers, Taylor-Green vortex, and lid-driven cavity. Notably, this method enables the discovery of stable, non-trivial LB parameter sets for the nonlinear case. Finally, to assess numerical stability or instability with these non-trivial parameters, a global stability analysis is conducted over the entire lattice, encompassing the boundary conditions of the linearized scheme.