The phenomenon of gravitational particle production can take place for quantum fields in curved spacetime. The abundance and energy spectrum of gravitationally produced particles is typically calculated by solving the field's mode equations on a timedependent background metric. For purposes of studying dark matter production in an inflationary cosmology, these mode equations are often solved numerically, which is computationally intensive, especially for the rapidly-oscillating high-momentum modes. However, these same modes are amenable to analytic evaluation via the Exact Wentzel-Kramers-Brillouin (EWKB) method, where gravitational particle production is a manifestation of the Stokes phenomenon. These analytic techniques have been used in the past to study gravitational particle production for spin-0 bosons. We extend the earlier work to study gravitational production of spin-1/2 and spin-3/2 fermions. We derive an analytic expression for the connection matrix (valid to all orders in perturbations) that relates Bogoliubov coefficients across a Stokes line connecting a merged pair of simple turning points. By comparing the analytic approximation with a direct numerical integration of the mode equations, we demonstrate an excellent agreement and highlight the utility of the Stokes phenomenon formalism applied to fermions. We discuss the implications for an analytic understanding of catastrophic particle production due to vanishing sound speed, which can occur for a spin-3/2 Rarita-Schwinger field.