Multifit

Horace (and the parent utilities library, Herbert) comes with a rich and powerful fitting syntax that is used to fit functions or models of S(Q,w) to one or more datasets. The documentation here is only meant to give an introduction and overview. For the full help, please use the Matlab documentation for the various fitting functions that can be obtained by using the doc command, for example doc sqw/multifit (fitting functions to sqw objects) or doc sqw/multifit_sqw (fitting S(Q,w) to sqw objects).

Multifit is one of the newer capabilities introduced into Horace - it was not included in the alpha test version of the software circulated until April 2009. The purpose of Multifit is to allow fitting of multiple data structures of any dimensionality to a global model. The information below is hopefully sufficiently comprehensive to allow you to work out how to use these functions, however you can also get help by looking at the Horace Masterclass file. If you are still having problems after consulting this, then get in touch with us at Horace Help.

Overview

There are two ways of fitting multiple dataset, one is to fit a global set of parameters, and the other is to fit each dataset individually. In other words the fits can be parallel or sequential. For the latter kind of fit you can use the functions fit_func and fit_sqw, with a first argument that is an array of dnd / sqw objects (all with the same dimensionality). In this case the fit program will run as if you had placed it inside a ‘’for’’ loop, and the resulting output fit object will be an array the same size as the input array.

However, in many cases what is required is to fit a global model to a dataset. For example, you may have taken a series of 2d slices that show some magnetic scattering, that were taken a special energies or wavevectors where non-magnetic scattering was minimal. Rather than fitting the slices one at a time and then working out some kind of average for the fit parameters, one can instead fit all of them to a single model with a single set of parameters. For this purpose you should use the command ’multifit’.

Additionally you may wish to keep your model function and background separate. For example, you may have several objects to fit for which the same model is appropriate, but for which the background varies significantly between models. In this case you can choose to fit the background for each object separately, or even use different background functions for different objects.

When fitting your model or background you may wish to keep the ratio of two fit parameters constant. This is possible for both the model function and the background function(s) when using multifit, by using what we have termed ‘’’binding’’’ of parameters.

There are three versions of this function. The first two are analogous to the regular fitting functions, one for full ‘’S(Q,$\omega$)’’ models and another for more general functions. These are called multifit_sqw and multifit_func respectively. In both cases the separate background functions are simply f(x,y,...) and are specific to the dimensionality of the data being fitted. There is additionally multifit_sqw_sqw, for which both the main function and the background function are given by ‘’S(Q,$\omega$)’’ models.

Finally note that the keywords that are used for the standard fit function also apply in exactly the same way in multifit.

Detail

In the discussion below we will start by considering the most basic syntax that can be provided to the function, and will then gradually introduce further complexity. These notes are essentially the same as the in-built help for multifit, albeit with a little more detail. To read this abridged version, type

>> help multifit_sqw


There are several different ways of calling the multifit function, with varying levels of complexity:

>> [wout, fitdata] = multifit_sqw (w, func, pin)
>> [wout, fitdata] = multifit_sqw (w, func, pin, pfree)
>> [wout, fitdata] = multifit_sqw (w, func, pin, pfree, pbind)


The inputs have a similar form to those used in fit. For the first case the array of objects w are fitted using the function func, with an initial guess of the parameters supplied in the vector or cell array pin. If pin is a vector then all of its elements are fittable parameters. If it is a cell array then its first cell must be a vector of this form, and the remaining cells contain additional information passed to the fit function, that cannot be varied (e.g. magnetic form factor lookup tables).

If no pfree is specified then it is assumed that all fit parameters can be varied. Otherwise pfree must be a vector with the same length as the number of fit parameters, and comprising 1s and 0s (for the parameter being free or fixed respectively).

pbind is an additional argument which is a halfway house between allowing parameters to be free or fixed. It allows you to specify that some parameters must maintain a fixed ratio, but that their absolute values can be varied. For example, you may wish to fix the ratio of two exchange constants, but allow their magnitude to vary in order to find the best fit.

There are several different things you can do with parameter binding, so you must be careful about how you write pbind. The simplest form for pbind is a 2-element cell array, such as pbind = {1,3}, which binds the first and third parameters to a constant ratio, determined by the ratio between the parameters given in pin. You can specify multiple pairs of bindings by saying pbind = {{1,3},{4,3},{5,6}}, so parameters 1 and 3, 4 and 3, and 5 and 6 are bound in the ratios of the initial guesses. You can also explicitly provide the ratios for the bindings (rather than using the ratio of the initial guess) by writing, for example, pbind = {{1,3,0,7.4},{4,3,0,0.023}} in which case parameter 1 and 3 are bound in the ratio 7.4, and 4 and 3 are bound with a ratio 0.023. Notice the extra zero – this is required so that the program knows what you are trying to do. Finally, you can choose to bind parameters in the fit function to parameters in the separate background function(s). The syntax for this would be, for example, pbind = {1,3,[1,2],3.1} - whereby parameter 1 is bound to background parameter 3 for the 2nd background function, and the ratio is 3.1.

>> [wout, fitdata] = multifit_sqw (..., bkdfunc, bpin)
>> [wout, fitdata] = multifit_sqw (..., bkdfunc, bpin, bpfree)
>> [wout, fitdata] = multifit_sqw (..., bkdfunc, bpin, bpfree, bpbind)


As mentioned above, when using multifit the background can be fitted separately. You can either specify a different background function for each object to be fitted, or use a single function. In each case the backgrounds are not fitted globally, so the background parameters for each fit can be different, even though the model fit is the same. The syntax for the background fitting is a little more complicated than for the main fitting, essentially because there are more permutations for keeping parameters fixed and/or bound.

There are two options for specifying the background function bkdfunc. You either specify a single function handle (in the same way as you did for the fitting function), or you provide a cell array of function handles. In the former case the same background function (but with different parameters) is used for every object to be fitted. In the latter case the cell array must have the same number of elements as there are objects to fit. Each element is then a function handle (e.g. @bg_func) which specifies the function that will be used to model the background for each object to be fitted. This means that even if you only wish to use a different background function for one of the objects to be fitted, you must give a cell array of function handles (for that case all of the function handles would be the same, except for one).

For specifying the initial guesses for the parameters passed to the background function(s) the syntax is also pretty similar to that used for the fitting function. If only a single background function is to be used then the parameters bpin can be specified as either a vector or a cell array. If the former then the elements of the vector are the parameters required by the background function. If the latter, then the first element of the cell array should be a vector of parameters, and the other elements are extra parameters (e.g. lookup tables, etc) that are required by the background function. If more than one background function is being used, or you wish to supply a different set of initial guesses for different objects to be fitted, then bpin must be a cell array with the same number of elements as there are objects to be fitted. The form of each element of the cell array would then be the same as that described above, i.e. a vector of parameters, or a cell array of parameters (so you could have a cell array of cell arrays…).

You can specify which parameters to vary or keep fixed for the background in a similar way to how you do it for the fitting function. You can either specify a vector of 1s and 0s if only one background function is being used. If more than one background function is used, or you do not wish to vary the parameters for every object, you can specify a cell array which has a vector of 1s and 0s in each element.

Finally for the most complicated part – specifying the bindings for the background. You should be especially careful here, as there could be cell arrays all over the place, and the precise format depends on a variety of factors. The way to lay this out is best illustrated by a series of examples of single binding descriptions:

1) bpbind = {1,4}. The background parameter bp(1) is bound to background parameter bp(4) with a fixed ratio, determined by the initial guess. This binding will apply for all objects and for all background functions. 2) bpbind = {2,3,[7,2]}. bp(2) is bound to bp(3) for background function bkdfunc{7,2} only, in the case where multiple functions have been specified. For all other background functions there is no binding. The ratio is determined from the initial guesses. 3) bpbind = {5,11,0}. bp(5) is bound to parameter p(11) of the global fitting function. This binding will apply for all objects and all background functions. The ratio is determined by the initial guess. 4) bpbind = {{1,4},{2,3,[7,2]},{5,11,0}}. Several different bindings (see above) are defined. 5) bpbind = {5,11,0,0.013}. Background parameter bp(5) is bound to global parameter p(11), with a ratio explicitly given as 0.013, overriding the initial guess. 6) bpbind = {1,4,[7,2],14.15}. bp(1) is bound to bp(4) for bg_func(7,2) only, with a ratio of 14.15, which overrides the ratio determined from the initial guess.

As an illustration of the need for care, consider the following:

bpbind1 = { {1,4}, {2,3,[7,2]}, {5,11,0} }
bpbind2 = { { {1,4}, {2,3,[7,2]}, {5,11,0} } }


bpbind1 specifies the bindings for 3 different backgrounds. bpbind2 specifies a set of bindings applicable to all backgrounds.

>> [wout, fitdata] = multifit_sqw (…,keyword, value)


Keywords relating to which parts of the data to fit, what information to display in the Matlab command window whilst the fit is running, etc, can be chosen. The keywords and values are the same as for the functions fit and fit_sqw.

The outputs from multifit are as follows. wout is an array of dnd / sqw objects that matches the input, but the intensity values are the results of the fitting. fitdata is a structure array with the following fields:

1) fitdata.p – the global parameter values. 2) fitdata.sig – the errors on the global parameters (=0 if a parameter was fixed). 3) fitdata.bp – the background parameters. This is either a vector, or a cell array of vectors, depending on how the input was arranged. 4) fitdata.bsig – the errors on the background parameters, in either a vector or a cell array of vectors. 5) fitdata.corr – the correlation matrix for the free parameters. 6) fitdata.chisq – the value of the reduced Chi^2 of the fit, i.e. divided by (number of data points minus number of free parameters). 7) fitdata.pnames – parameter names (if they were specified in the fit function). 8) fitdata.bpnames – background parameter names (if they were specified in the background function(s)).

Examples

Let us illustrate how multifit is actually used by means of a series of examples. For all of the following we shall assume that we wish to fit 3 two-dimensional sqw objects comprising constant energy slices in the (h,0,0),(0,k,0) plane, and that we wish to fit a quartet of Gaussian profile peaks, much as we did in the Getting Started section of this manual.

No background functions

In this case we allow all of the parameters to vary, and we bind the second and third parameters to a ratio of 1.

win(1)=first_obj;
win(2)=second_obj;
win(3)=third_obj;

pin=[6 1 1 0.1 1.25 6];

pfree=ones(size(pin{1}));

pbind={2,3,0,1};

[wout, fitdata] = multifit_sqw (win, @demo_4gauss, pin, pfree, pbind);


Simple background functions

In this case we will have the same as above, but will add a single background function that we allow to vary between the objects to be fitted. We supply different initial guesses for the background parameters, and we bind the third background parameter to the 1st background parameter in the ratio specified by the initial guess.

pin=[6 1 1 0.1 1.25 6];

pfree=ones(size(pin{1}));

pbind={2,3,0,1};

bpin=cell(1,3);
bpin{1}=[1 1 0];
bpin{2}=[1 1.5 0.2];
bpin{3}=[2 -0.1 1];

bpfree=[1,1,1];

bpbind={3,1};

[wout, fitdata] = multifit_sqw (win, @demo_4gauss, pin, pfree, pbind, @sloping_background,...
bpin, bpfree, bpbind);


The Fully Monty

Let us now demonstrate a more complicated use of multifit. In this case we will have a different background function for each object. The third of these background functions will require additional information passed to it (a lookup table), necessitating the use of a nested cell array. We will bind the third parameter to the first parameter of the first background function, and likewise for the second background function. Finally we will also bind the second parameter of the background functions to the first parameter of the fitting function.

pin=[6 1 1 0.1 1.25 6];

pfree=ones(size(pin{1}));

pbind={2,3,0,1};

bg_funcs=cell(1,3);
bg_funcs{1}=@constant_background;
bg_funcs{2}=@sloping_background;
bg_funcs{3}=@special_background;

bpin=cell(1,3);
bpin{1}=[1 1 0];
bpin{2}=[1 1.5 0.2];
bg_cell=cell(1,2)
bg_cell{1}=[2 -0.1 1];
bg_cell{2}='C:\mprogs\lookup_table.dat';
bpin{3}=bg_cell;

bpfree=cell{1,3};
bpfree{1}=[1,1,1];
bpfree{2}=[1,1,0];
bpfree{3}=[0,1,0];

bpbind={{3,1,[1,1]},{3,1,[1,2]},{2,1,0}};

[wout, fitdata] = multifit_sqw (win, @demo_4gauss, pin, pfree, pbind, bg_funcs,...
bpin, bpfree, bpbind);


Semi-global fits

Horace allows you to fit several objects independently (i.e. sequentially) using the fit_sqw and fit functions. You can also fit several objects with a global set of parameters using multifit. However there is a third thing you might wish to do, halfway between these two, which is to fit a model where some of the parameters are global but others are object-specific. For example, you may wish to fit several 1d cuts which show some magnetic scattering, but are taken in different Brillouin zones so that the intensity scaling is very different. If the magnetic form-factor alone does not describe the intensity differences you may wish to have a different intensity scale for each cut, but to keep other parameters global. Previously there is no straightforward way of doing this in multifit, however it is now possible to do this using the function multifit_sqw_sqw.

multifit_sqw_sqw allows you to have sqw models for both the spectral function and the background function. You can then bind parameters in the background function to parameters in the spectral function, i.e. the background function becomes the spectral function.

In order to make some parameters global, and others object-specific, you can be clever with the background bindings. Let us illustrate this with an example…

Suppose the spectral weight can be specified by a model which has I0, k0 and G as its input parameters, and that the background is just a constant B. We wish to allow I0 and B to be different for different objects, but k0 and G to be global. We specify a spectral function with the handle @null_func, the Matlab code for which would be something like:

function weight = null_func(qh,qk,ql,en,p)

weight=zeros(size(qh));


We then write a Matlab function called real_func which takes I0, k0, G and B as input parameters and calculates the spectral weight plus background. We then use multifit as follows:

[wfit,fitdata] = multifit_sqw_sqw( win,@null_func,[k0_guess,G_guess],[1,1],…
@real_func,[I0_guess,k0_guess,G_guess,B_guess], [1,1,1,1],…
{  { {2,1,0,1}, {3,2,0,1} }  } );


The binding argument we have specified applies to all backgrounds. We have bound parameter 2 of real_func to parameter 1 of null_func, in the ratio 1:1; and we have bound parameter 3 of real_func to parameter 2 of null_func, in the ratio 1:1.

When this code is executed the actual fitting all happens with real_func, since null_func does not do anything, however by passing k0 and G to null_func we can ensure that they are global parameters. Of course if your model is complicated then writing the function real_func might be rather involved, but it can be done if you have the patience!