Graph signal processing (GSP) advances spectral analysis on irregular domains. However, existing two-dimensional graph fractional Fourier transform (2D-GFRFT) employs a single fractional order for both factor graphs, thereby limiting its adaptability to heterogeneous signals. We proposed the two-dimensional graph bi-fractional Fourier transform (2D-GBFRFT), which assigns independent fractional orders to the factor graphs of a Cartesian product while preserving separability. We established invertibility, unitarity, and index additivity, and developed two filtering schemes: a Wiener-style design through grid search and a differentiable framework that jointly optimizes transform orders and diagonal spectral filters. We further introduced a hybrid interpolation with the joint time-vertex fractional Fourier transform (JFRFT), controlled by a tunable parameter that balances the two methods. In the domains of synthetic Cartesian product graph signals, authentic temporal graph datasets, and dynamic image deblurring, 2D-GBFRFT consistently surpasses 2D-GFRFT and enhances JFRFT. Experimental results confirmed the versatility and superior performance of 2D-GBFRFT for filtering in GSP.