  Vagus Nerve Stimulation (VNS) is an established palliative treatment for drug
resistant epilepsy. While effective for many patients, its mechanism of action
is incompletely understood. Predicting individuals' response, or optimum
stimulation parameters, is challenging. Computational modelling has informed
other problems in epilepsy but, to our knowledge, has not been applied to VNS.
  We started with an established, four-population neural mass model (NMM),
capable of reproducing the seizure-like dynamics of a thalamocortical circuit.
We extended this to include 18 further neural populations, representing nine
other brain regions relevant to VNS, with connectivity based on existing
literature. We modelled stimulated afferent vagal fibres as projecting to the
nucleus tractus solitarius (NTS), which receives input from the vagus nerve in
vivo.
  Bifurcation analysis of a deterministic version of the model showed higher
background NTS input made the model monostable at a fixed point (FP),
representing normal activity, while lower inputs produce bistability between
the FP and a limit cycle (LC), representing the seizure state.
  Adding noise produced transitions between seizure and normal states. This
stochastic model spent decreasing time in the seizure state with increasing
background NTS input, until seizures were abolished, consistent with the
deterministic model.
  Simulated VNS stimulation, modelled as a 30 Hz square wave, was summed with
the background input to the NTS and was found to reduce total seizure duration
in a dose-dependent manner, similar to expectations in vivo.
  We have successfully produced an in silico model of VNS in epilepsy,
capturing behaviour seen in vivo. This may aid understanding therapeutic
mechanisms of VNS in epilepsy and provides a starting point to (i) determine
which patients might respond best to VNS, and (ii) optimise individuals'
treatments.
