library: libMinuit #include "TLinearFitter.h"

# TLinearFitter

class description - source file - inheritance tree (.pdf)

## class TLinearFitter : public TVirtualFitter

| Inheritance Chart:
`TObject`
| <- | `TNamed`
| <- | `TVirtualFitter`
| <- | `TLinearFitter` |

        private:

              void AddToDesign(Double_t* x, Double_t y, Double_t e)
              void CreateSubset(Int_t ntotal, Int_t h, Int_t* index)
          Double_t CStep(Int_t step, Int_t h, Double_t* residuals, Int_t* index, Int_t* subdat, Int_t start, Int_t end)
              void Graph2DLinearFitter(Double_t h)
              void GraphLinearFitter(Double_t h)
              void HistLinearFitter()
          Double_t KOrdStat(Int_t ntotal, Double_t* a, Int_t k, Int_t* work)
            Bool_t Linf()
              void MultiGraphLinearFitter(Double_t h)
             Int_t Partition(Int_t nmini, Int_t* indsubdat)
              void RDraw(Int_t* subdat, Int_t* indsubdat)

        public:

                            TLinearFitter()
                            TLinearFitter(Int_t ndim, const char* formula, Option_t* opt = "D")
                            TLinearFitter(Int_t ndim)
                            TLinearFitter(TFormula* function, Option_t* opt = "D")
                    virtual ~TLinearFitter()
               virtual void AddPoint(Double_t* x, Double_t y, Double_t e = 1)
               virtual void AssignData(Int_t npoints, Int_t xncols, Double_t* x, Double_t* y, Double_t* e = 0)
               virtual void Chisquare()
           virtual Double_t Chisquare(Int_t, Double_t*) const
             static TClass* Class()
               virtual void Clear(Option_t* option = "")
               virtual void ClearPoints()
               virtual void Eval()
               virtual void EvalRobust(Double_t h = -1)
              virtual Int_t ExecuteCommand(const char* command, Double_t* args, Int_t nargs)
               virtual void FixParameter(Int_t ipar)
               virtual void FixParameter(Int_t ipar, Double_t parvalue)
           virtual Double_t GetChisquare()
          virtual Double_t* GetCovarianceMatrix() const
               virtual void GetCovarianceMatrix(TMatrixD& matr)
           virtual Double_t GetCovarianceMatrixElement(Int_t i, Int_t j) const
               virtual void GetErrors(TVectorD& vpar)
              virtual Int_t GetErrors(Int_t, Double_t&, Double_t&, Double_t&, Double_t&) const
               virtual void GetFitSample(TBits& bits)
              virtual Int_t GetNumberFreeParameters() const
              virtual Int_t GetNumberTotalParameters() const
           virtual Double_t GetParameter(Int_t ipar) const
              virtual Int_t GetParameter(Int_t, char*, Double_t&, Double_t&, Double_t&, Double_t&) const
               virtual void GetParameters(TVectorD& vpar)
           virtual Double_t GetParError(Int_t ipar) const
           virtual Double_t GetParSignificance(Int_t ipar) const
           virtual Double_t GetParTValue(Int_t ipar) const
              virtual Int_t GetStats(Double_t&, Double_t&, Double_t&, Int_t&, Int_t&) const
           virtual Double_t GetSumLog(Int_t)
            virtual TClass* IsA() const
             virtual Bool_t IsFixed(Int_t ipar) const
               virtual void PrintResults(Int_t level, Double_t amin = 0) const
               virtual void ReleaseParameter(Int_t ipar)
               virtual void SetDim(Int_t n)
               virtual void SetFitMethod(const char*)
               virtual void SetFormula(const char* formula)
               virtual void SetFormula(TFormula* function)
              virtual Int_t SetParameter(Int_t, const char*, Double_t, Double_t, Double_t, Double_t)
               virtual void ShowMembers(TMemberInspector& insp, char* parent)
               virtual void StoreData(Bool_t store)
               virtual void Streamer(TBuffer& b)
                       void StreamerNVirtual(TBuffer& b)
             virtual Bool_t UpdateMatrix()

### Data Members

        private:

             TVectorD fParams         vector of parameters
          TMatrixDSym fParCovar       matrix of parameters' covariances
             TVectorD fTValues        T-Values of parameters
             TVectorD fParSign        significance levels of parameters
          TMatrixDSym fDesign         matrix AtA
          TMatrixDSym fDesignTemp     ! temporary matrix, used for num.stability
          TMatrixDSym fDesignTemp2    !
          TMatrixDSym fDesignTemp3    !
             TVectorD fAtb            vector Atb
             TVectorD fAtbTemp        ! temporary vector, used for num.stability
             TVectorD fAtbTemp2       !
             TVectorD fAtbTemp3       !
              Bool_t* fFixedParams    array of fixed/released params
            TObjArray fFunctions      array of basis functions
             TVectorD fY              the values being fit
             Double_t fY2             sum of square of y, used for chisquare
             Double_t fY2Temp         ! temporary variable used for num.stability
             TMatrixD fX              values of x
             TVectorD fE              the errors if they are known
            TFormula* fInputFunction  the function being fit
                Int_t fNpoints        number of points
                Int_t fNfunctions     number of basis functions
                Int_t fFormulaSize    length of the formula
                Int_t fNdim           number of dimensions in the formula
                Int_t fNfixed         number of fixed parameters
                Int_t fSpecial        =100+n if fitting a polynomial of deg.n
                char* fFormula        the formula
               Bool_t fIsSet          Has the formula been set?
               Bool_t fStoreData      Is the data stored?
             Double_t fChisquare      Chisquare of the fit
                Int_t fH              number of good points in robust fit
               Bool_t fRobust         true when performing a robust fit
                TBits fFitsample      indices of points, used in the robust fit

## Class Description

     The Linear Fitter - fitting functions that are LINEAR IN PARAMETERS

     Linear fitter is used to fit a set of data points with a linear
     combination of specified functions. Note, that "linear" in the name
     stands only for the model dependency on parameters, the specified
     functions can be nonlinear.
     The general form of this kind of model is

              y(x) = a[0] + a[1]*f[1](x)+...a[n]*f[n](x)

     Functions f are fixed functions of x. For example, fitting with a
     polynomial is linear fitting in this sense.

                             The fitting method

     The fit is performed using the Normal Equations method with Cholesky
     decomposition.

                             Why should it be used?

     The linear fitter is considerably faster than general non-linear
     fitters and doesn't require to set the initial values of parameters.

                              Using the fitter:

     1.Adding the data points:
      1.1 To store or not to store the input data?
          - There are 2 options in the constructor - to store or not
            store the input data. The advantages of storing the data
            are that you'll be able to reset the fitting model without
            adding all the points again, and that for very large sets
            of points the chisquare is calculated more precisely.
            The obvious disadvantage is the amount of memory used to
            keep all the points.
          - Before you start adding the points, you can change the
            store/not store option by StoreData() method.
      1.2 The data can be added:
          - simply point by point - AddPoint() method
          - an array of points at once:
            If the data is already stored in some arrays, this data
            can be assigned to the linear fitter without physically
            coping bytes, thanks to the Use() method of
            TVector and TMatrix classes - AssignData() method

     2.Setting the formula
      2.1 The linear formula syntax:
          -Additive parts are separated by 2 plus signes "++"
           --for example "1 ++ x" - for fitting a straight line
          -All standard functions, undrestood by TFormula, can be used
           as additive parts
           --TMath functions can be used too
          -Functions, used as additive parts, shouldn't have any parameters,
           even if those parameters are set.
           --for example, if normalizing a sum of a gaus(0, 1) and a
             gaus(0, 2), don't use the built-in "gaus" of TFormula,
             because it has parameters, take TMath::Gaus(x, 0, 1) instead.
          -Polynomials can be used like "pol3", .."polN"
          -If fitting a more than 3-dimensional formula, variables should
           be numbered as follows:
           -- x0, x1, x2... For example, to fit  "1 ++ x0 ++ x1 ++ x2 ++ x3*x3"
      2.2 Setting the formula:
        2.2.1 If fitting a 1-2-3-dimensional formula, one can create a
              TF123 based on a linear expression and pass this function
              to the fitter:
              --Example:
                TLinearFitter *lf = new TLinearFitter();
                TF2 *f2 = new TF2("f2", "x ++ y ++ x*x*y*y", -2, 2, -2, 2);
                lf->SetFormula(f2);
              --The results of the fit are then stored in the function,
                just like when the TH1::Fit or TGraph::Fit is used
              --A linear function of this kind is by no means different
                from any other function, it can be drawn, evaluated, etc.
        2.2.2 There is no need to create the function if you don't want to,
              the formula can be set by expression:
              --Example:
                // 2 is the number of dimensions
                TLinearFitter *lf = new TLinearFitter(2);
                lf->SetFormula("x ++ y ++ x*x*y*y");
              --That's the only way to go, if you want to fit in more
                than 3 dimensions
        2.2.3 The fastest functions to compute are polynomials and hyperplanes.
              --Polynomials are set the usual way: "pol1", "pol2",...
              --Hyperplanes are set by expression "hyp3", "hyp4", ...
              ---The "hypN" expressions only work when the linear fitter
                 is used directly, not through TH1::Fit or TGraph::Fit.
                 To fit a graph or a histogram with a hyperplane, define
                 the function as "1++x++y".
              ---A constant term is assumed for a hyperplane, when using
                 the "hypN" expression, so "hyp3" is in fact fitting with
                 "1++x++y++z" function.
              --Fitting hyperplanes is much faster than fitting other
                expressions so if performance is vital, calculate the
                function values beforehand and give them to the fitter
                as variables
              --Example:
                You want to fit "sin(x)|cos(2*x)" very fast. Calculate
                sin(x) and cos(2*x) beforehand and store them in array *data.
                Then:
                TLinearFitter *lf=new TLinearFitter(2, "hyp2");
                lf->AssignData(npoint, 2, data, y);

      2.3 Resetting the formula
        2.3.1 If the input data is stored (or added via AssignData() function),
              the fitting formula can be reset without re-adding all the points.
              --Example:
                TLinearFitter *lf=new TLinearFitter("1++x++x*x");
                lf->AssignData(n, 1, x, y, e);
                lf->Eval()
                //looking at the parameter significance, you see,
                // that maybe the fit will improve, if you take out
                // the constant term
                lf->SetFormula("x++x*x");
                lf->Eval();
                ...
        2.3.2 If the input data is not stored, the fitter will have to be
              cleared and the data will have to be added again to try a
              different formula.

     3.Accessing the fit results
      3.1 There are methods in the fitter to access all relevant information:
          --GetParameters, GetCovarianceMatrix, etc
          --the t-values of parameters and their significance can be reached by
            GetParTValue() and GetParSignificance() methods
      3.2 If fitting with a pre-defined TF123, the fit results are also
          written into this function.

     4.Robust fitting - Least Trimmed Squares regression (LTS)
       Outliers are atypical(by definition), infrequant observations; data points
       which do not appear to follow the characteristic distribution of the rest
       of the data. These may reflect genuine properties of the underlying
       phenomenon(variable), or be due to measurement errors or anomalies which
       shouldn't be modelled. (StatSoft electronic textbook)

       Even a single gross outlier can greatly influence the results of least-
       squares fitting procedure, and in this case use of robust(resistant) methods
       is recommended.

       The method implemented here is based on the article and algorithm:
       "Computing LTS Regression for Large Data Sets" by
       P.J.Rousseeuw and Katrien Van Driessen
       The idea of the method is to find the fitting coefficients for a subset
       of h observations (out of n) with the smallest sum of squared residuals.
       The size of the subset h should lie between (npoints + nparameters +1)/2
       and n, and represents the minimal number of good points in the dataset.
       The default value is set to (npoints + nparameters +1)/2, but of course
       if you are sure that the data contains less outliers it's better to change
       h according to your data.

       To perform a robust fit, call EvalRobust() function instead of Eval() after
       adding the points and setting the fitting function.
       Note, that standard errors on parameters are not computed!

TLinearFitter()

    default c-tor, input data is stored
    If you don't want to store the input data,
    run the function StoreData(kFALSE) after constructor

TLinearFitter(Int_t ndim)

    The parameter stands for number of dimensions in the fitting formula
    The input data is stored. If you don't want to store the input data,
    run the function StoreData(kFALSE) after constructor

TLinearFitter(Int_t ndim, const char *formula, Option_t *opt)

    First parameter stands for number of dimensions in the fitting formula
    Second parameter is the fitting formula: see class description for formula syntax
    Options:
    The option is to store or not to store the data
    If you don't want to store the data, choose "" for the option, or run
    StoreData(kFalse) member function after the constructor

TLinearFitter(TFormula *function, Option_t *opt)

    This constructor uses a linear function. How to create it?
    TFormula now accepts formulas of the following kind:
    TFormula("f", "x++y++z++x*x"). Other than the look, it's in no
    way different from the regular formula, it can be evaluated,
    drawn, etc.
    The option is to store or not to store the data
    If you don't want to store the data, choose "" for the option, or run
    StoreData(kFalse) member function after the constructor

~TLinearFitter()

     Linear fitter cleanup.

void AddPoint(Double_t *x, Double_t y, Double_t e)

    Adds 1 point to the fitter.
    First parameter stands for the coordinates of the point, where the function is measured
    Second parameter - the value being fitted
    Third parameter - weight(measurement error) of this point (=1 by default)

void AssignData(Int_t npoints, Int_t xncols, Double_t *x, Double_t *y, Double_t *e)

    This function is to use when you already have all the data in arrays
    and don't want to copy them into the fitter. In this function, the Use() method
    of TVectorD and TMatrixD is used, so no bytes are physically moved around.
    First parameter - number of points to fit
    Second parameter - number of variables in the model
    Third parameter - the variables of the model, stored in the following way:
    (x0(0), x1(0), x2(0), x3(0), x0(1), x1(1), x2(1), x3(1),...

void AddToDesign(Double_t *x, Double_t y, Double_t e)

    Add a point to the AtA matrix and to the Atb vector.

void Clear(Option_t * /*option*/)

    Clears everything. Used in TH1::Fit and TGraph::Fit().

void ClearPoints()

    To be used when different sets of points are fitted with the same formula.

void Chisquare()

    Calculates the chisquare.

void Eval()

     Perform the fit and evaluate the parameters

void FixParameter(Int_t ipar)

    Fixes paramter #ipar at its current value.

void FixParameter(Int_t ipar, Double_t parvalue)

    Fixes parameter #ipar at value parvalue.

void ReleaseParameter(Int_t ipar)

    Releases parameter #ipar.

Double_t GetChisquare()

     Get the Chisquare.

void GetCovarianceMatrix(TMatrixD &matr)

void GetErrors(TVectorD &vpar)

void GetParameters(TVectorD &vpar)

Double_t GetParError(Int_t ipar) const

void GetFitSample(TBits &bits)

void SetDim(Int_t ndim)

    set the number of dimensions

void SetFormula(const char *formula)

    Additive parts should be separated by "++".
    Examples (ai are parameters to fit):
    1.fitting function: a0*x0 + a1*x1 + a2*x2
      input formula "x0++x1++x2"
    2.TMath functions can be used:
      fitting function: a0*TMath::Gaus(x0, 0, 1) + a1*x1
      input formula:    "TMath::Gaus(x0, 0, 1)++x1"
    fills the array of functions

void SetFormula(TFormula *function)

    Set the fitting function.

Bool_t UpdateMatrix()

Int_t ExecuteCommand(const char *command, Double_t *args, Int_t /*nargs*/)

    To use in TGraph::Fit and TH1::Fit().

void PrintResults(Int_t level, Double_t /*amin*/) const

     Level = 3 (to be consistent with minuit)  prints parameters and parameter
     errors.

void GraphLinearFitter(Double_t h)

    Used in TGraph::Fit().

void Graph2DLinearFitter(Double_t h)

void MultiGraphLinearFitter(Double_t h)

void HistLinearFitter()

     Minimization function for H1s using a Chisquare method.

void EvalRobust(Double_t h)

    Finds the parameters of the fitted function in case data contains
    outliers.
    Parameter h stands for the minimal fraction of good points in the
    dataset (h < 1, i.e. for 70% of good points take h=0.7).
    The default value of h*Npoints is  (Npoints + Nparameters+1)/2
    If the user provides a value of h smaller than above, default is taken
    See class description for the algorithm details

void CreateSubset(Int_t ntotal, Int_t h, Int_t *index)

    Creates a p-subset to start
    ntotal - total number of points from which the subset is chosen

Double_t CStep(Int_t step, Int_t h, Double_t *residuals, Int_t *index, Int_t *subdat, Int_t start, Int_t end)

    The CStep procedure, as described in the article

Bool_t Linf()

Int_t Partition(Int_t nmini, Int_t *indsubdat)

    divides the elements into approximately equal subgroups
    number of elements in each subgroup is stored in indsubdat
    number of subgroups is returned

void RDraw(Int_t *subdat, Int_t *indsubdat)

    Draws ngroup nonoverlapping subdatasets out of a dataset of size n
    such that the selected case numbers are uniformly distributed from 1 to n

Double_t KOrdStat(Int_t ntotal, Double_t *a, Int_t k, Int_t *work)

    copy of the TMath::KOrdStat because I need an Int_t work array

### Inline Functions

                   void GetCovarianceMatrix(TMatrixD& matr)
               Double_t GetCovarianceMatrixElement(Int_t i, Int_t j) const
                  Int_t GetNumberTotalParameters() const
                  Int_t GetNumberFreeParameters() const
               Double_t GetParameter(Int_t ipar) const
               Double_t GetParTValue(Int_t ipar) const
               Double_t GetParSignificance(Int_t ipar) const
                 Bool_t IsFixed(Int_t ipar) const
                   void StoreData(Bool_t store)
               Double_t Chisquare(Int_t, Double_t*) const
                  Int_t GetErrors(Int_t, Double_t&, Double_t&, Double_t&, Double_t&) const
                  Int_t GetParameter(Int_t, char*, Double_t&, Double_t&, Double_t&, Double_t&) const
                  Int_t GetStats(Double_t&, Double_t&, Double_t&, Int_t&, Int_t&) const
               Double_t GetSumLog(Int_t)
                   void SetFitMethod(const char*)
                  Int_t SetParameter(Int_t, const char*, Double_t, Double_t, Double_t, Double_t)
                TClass* Class()
                TClass* IsA() const
                   void ShowMembers(TMemberInspector& insp, char* parent)
                   void Streamer(TBuffer& b)
                   void StreamerNVirtual(TBuffer& b)

Author: Anna Kreshuk 04/03/2005
Last update: root/minuit:$Name: $:$Id: TLinearFitter.cxx,v 1.11 2005/06/23 10:04:08 brun Exp $
Copyright (C) 1995-2005, Rene Brun and Fons Rademakers. *