Bilevel optimization involves a hierarchy of two optimization problems. Some variables in the upper level problem are constrained to be optimal in the lower level problem. We investigate a connection between bilevel optimization and stability of nonlinear dynamical systems. An important system characteristic is the size of its region of attraction. One can estimate the region of attraction through the sublevel sets of a Lyapunov function. We argue that this problem in its general form is indeed a bilevel optimization problem. As a result, we arrive at two conclusions. First, our problem in its general form is NP-hard, since bilevel programming has been proved to be NP-hard even in its simplest form. Second, one can use algorithms developed for bilevel programming to estimate the region of attraction, if they can handle nonconvexity of the lower level problem.