Logo ROOT  
Reference Guide
fitLinearRobust.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_fit
3/// \notebook -js
4/// This tutorial shows how the least trimmed squares regression,
5/// included in the TLinearFitter class, can be used for fitting
6/// in cases when the data contains outliers.
7/// Here the fitting is done via the TGraph::Fit function with option "rob":
8/// If you want to use the linear fitter directly for computing
9/// the robust fitting coefficients, just use the TLinearFitter::EvalRobust
10/// function instead of TLinearFitter::Eval
11///
12/// \macro_image
13/// \macro_output
14/// \macro_code
15///
16/// \author Anna Kreshuk
17
18#include "TRandom.h"
19#include "TGraphErrors.h"
20#include "TF1.h"
21#include "TCanvas.h"
22#include "TLegend.h"
23
24void fitLinearRobust()
25{
26 //First generate a dataset, where 20% of points are spoiled by large
27 //errors
28 Int_t npoints = 250;
29 Int_t fraction = Int_t(0.8*npoints);
30 Double_t *x = new Double_t[npoints];
31 Double_t *y = new Double_t[npoints];
32 Double_t *e = new Double_t[npoints];
33 TRandom r;
34 Int_t i;
35 for (i=0; i<fraction; i++){
36 //the good part of the sample
37 x[i]=r.Uniform(-2, 2);
38 e[i]=1;
39 y[i]=1 + 2*x[i] + 3*x[i]*x[i] + 4*x[i]*x[i]*x[i] + e[i]*r.Gaus();
40 }
41 for (i=fraction; i<npoints; i++){
42 //the bad part of the sample
43 x[i]=r.Uniform(-1, 1);
44 e[i]=1;
45 y[i] = 1 + 2*x[i] + 3*x[i]*x[i] + 4*x[i]*x[i]*x[i] + r.Landau(10, 5);
46 }
47
48 TGraphErrors *grr = new TGraphErrors(npoints, x, y, 0, e);
49 grr->SetMinimum(-30);
50 grr->SetMaximum(80);
51 TF1 *ffit1 = new TF1("ffit1", "pol3", -5, 5);
52 TF1 *ffit2 = new TF1("ffit2", "pol3", -5, 5);
53 ffit1->SetLineColor(kBlue);
54 ffit2->SetLineColor(kRed);
55 TCanvas *myc = new TCanvas("myc", "Linear and robust linear fitting");
56 myc->SetGrid();
57 grr->Draw("ap");
58 //first, let's try to see the result sof ordinary least-squares fit:
59 printf("Ordinary least squares:\n");
60 grr->Fit(ffit1);
61 //the fitted function doesn't really follow the pattern of the data
62 //and the coefficients are far from the real ones
63
64 printf("Resistant Least trimmed squares fit:\n");
65 //Now let's try the resistant regression
66 //The option "rob=0.75" means that we want to use robust fitting and
67 //we know that at least 75% of data is good points (at least 50% of points
68 //should be good to use this algorithm). If you don't specify any number
69 //and just use "rob" for the option, default value of (npoints+nparameters+1)/2
70 //will be taken
71 grr->Fit(ffit2, "+rob=0.75");
72 //
73 TLegend *leg = new TLegend(0.6, 0.8, 0.89, 0.89);
74 leg->AddEntry(ffit1, "Ordinary least squares", "l");
75 leg->AddEntry(ffit2, "LTS regression", "l");
76 leg->Draw();
77
78 delete [] x;
79 delete [] y;
80 delete [] e;
81
82}
ROOT::R::TRInterface & r
Definition: Object.C:4
#define e(i)
Definition: RSha256.hxx:103
int Int_t
Definition: RtypesCore.h:43
double Double_t
Definition: RtypesCore.h:57
@ kRed
Definition: Rtypes.h:64
@ kBlue
Definition: Rtypes.h:64
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition: TAttLine.h:40
The Canvas class.
Definition: TCanvas.h:27
1-Dim function class
Definition: TF1.h:210
A TGraphErrors is a TGraph with error bars.
Definition: TGraphErrors.h:26
virtual void SetMaximum(Double_t maximum=-1111)
Set the maximum of the graph.
Definition: TGraph.cxx:2251
virtual TFitResultPtr Fit(const char *formula, Option_t *option="", Option_t *goption="", Axis_t xmin=0, Axis_t xmax=0)
Fit this graph with function with name fname.
Definition: TGraph.cxx:1064
virtual void Draw(Option_t *chopt="")
Draw this graph with its current attributes.
Definition: TGraph.cxx:760
virtual void SetMinimum(Double_t minimum=-1111)
Set the minimum of the graph.
Definition: TGraph.cxx:2260
This class displays a legend box (TPaveText) containing several legend entries.
Definition: TLegend.h:23
virtual void SetGrid(Int_t valuex=1, Int_t valuey=1)
Definition: TPad.h:330
This is the base class for the ROOT Random number generators.
Definition: TRandom.h:27
Double_t y[n]
Definition: legend1.C:17
Double_t x[n]
Definition: legend1.C:17
leg
Definition: legend1.C:34