Designing efficient and user-friendly multimodal transit networks is critical for modern urban mobility. We study a novel stochastic multimodal transit network design problem that integrates fixed-route services with on-demand shuttles, explicitly accounting for heterogeneous rider preferences, uncertain travel times, and passenger demand. The hierarchical decision-making process is modeled using a two-stage stochastic bilevel optimization problem, where the transit agency (leader) determines the network design, and riders (followers) select their preferred routes based on realized traffic conditions. The model inherits the complexity of a nonconvex bilevel structure with stochastic programming, posing significant computational challenges.
To address this, we first develop an equivalent single-level mixed integer linear programming (MILP) reformulation by introducing a response search algorithm that efficiently enumerates critical follower route choices. To further enhance scalability, we propose a decomposition method that combines a relaxed formulation with a subset of follower responses and iteratively strengthens it with valid cutting planes. Computational experiments on instances derived from a public transit network in Dalian, China, demonstrate the efficiency and effectiveness of our approaches, achieving over 10 times speedups compared to existing single-level reformulations. Additionally, a comprehensive case study on the Ann Arbor/Ypsilanti region in Michigan highlights practical benefits of the proposed framework, yielding up to 12% cost savings and up to 7% improvements in route convenience, demonstrating the value of the proposed stochastic bilevel model over deterministic or single-level counterparts.