How can I define my own fit statistic in Sherpa?

The load_user_stat function provides an interface for adding user-defined statistic functions in Sherpa. It accommodates user-defined functions for a statistic and statistical errors, in addition to defining a list of model parameters and hyperparameters for prior distributions (if prior desired).

For example, a simple user-defined statistic in Python Sherpa would look like this:

load_user_stat(statname, calc_stat_func, calc_err_func=None, priors={})

sherpa> def my_stat_func(data, model, staterror, syserror=None, weight=None, bkg=None): 
    # A simple function to replicate χ2
    fvec = ((data - model) / staterror)**2
    stat = fvec.sum()
    return (stat, fvec)

sherpa> def my_staterr_func(data):
    # A simple staterror function
    return numpy.sqrt(data)

sherpa> load_user_stat("mystat", my_stat_func, my_staterr_func)

sherpa> set_stat(mystat)

And a more complex user-defined statistic, with prior distributions, would look like this:

sherpa> def my_stat_func(data, model, staterror=None, syserror=None, weight=None, bkg=None):
    ...
    return (stat, fvec)

sherpa> load_user_stat("mystat", my_stat_func, priors=dict(mugamma=0.017, ..., gamma=abs1.nh, ... ))

sherpa> set_stat(mystat)

In these examples, 'stat' is a scalar statistic value and 'fvec' is an array of the statistic contributions of each bin. This method caters to both types of optimization methods in Sherpa: levmar expects the statistic contribution per bin, whereas simplex and moncar expect a single scalar statistic value.