  This paper extends the existing fractional Hawkes process to better model
mainshock-aftershock sequences of earthquakes. The fractional Hawkes process is
a self-exciting point process model with temporal decay kernel being a
Mittag-Leffler function. A maximum likelihood estimation scheme is developed
and its consistency is checked. It is then compared to the ETAS model on three
earthquake sequences in Southern California. The fractional Hawkes process
performs favourably against the ETAS model. Additionally, two parameters in the
fractional Hawkes process may have a fixed geophysical meaning dependent on the
study zone and the stage of the seismic cycle the zone is in.
