Why don't we use Fermi-Dirac statistics for deriving Richardson equation for thermionic emission? Why does the Richardson-Dushman equation need Maxwell-Boltzmann statistics instead of Fermi-Dirac statistics to model thermionic emission, though only fermions are our only interest in those emissions?
 A: In the original derivation published by Richardson in 1901, he assumed that only electrons in the thermal or Maxwell–Boltzmann tail of the distribution function $f(E)$ contribute to the current density. In other words, the electrons with higher energies $(E \gg E_F)$ can escape from metal. 
Fermi-Dirac statistics approaches Maxwell-Boltzmann statistics for high temperature, low densities and also when $E \gg E_f$ often referred to as the classical limit.
$$f(E) = \dfrac{1}{1+\exp\left(\dfrac{E-E_F}{kT}\right)} \approx \exp\left(-\frac{E-E_F}{kT}\right)$$
The efforts to understand statistical methods and the MB distribution enable an understanding of first formulation of thermionic emission, in which electrons were modeled as evaporating off the surface of a metal. Richardson’s original derivation applied a Carnot cycle to the electron gas in equilibrium with a metal. Dushman reconsidered Richardson’s equation and his expression for the coefficient agrees with more modern treatments apart from a factor of $2$ associated with electron spin (a purely quantum mechanical feature of the electron revealed by the Stern–Gerlach experiment in 1922). It is a pleasingly different method of determining current density from thermionic emission using thermodynamic arguments without invoking the quantum mechanical ideas that form the backbone of most modern treatments.
