The reason Stan, Nimble, PyMC fit via simulation is because joint distributions generally don’t have an analytical form and require global fitting methods to identify.<p>I’d hazard to guess that if this package fits incrementally without simulation that either it doesn’t aim for a global optimum or that it places substantial limits on the form models can take.