Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TGraphQQ.cxx
Go to the documentation of this file.
1// @(#)root/graf:$Id$
2// Author: Anna Kreshuk 18/11/2005
3
4/*************************************************************************
5 * Copyright (C) 1995-2005, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12#include "TGraphQQ.h"
13#include "TAxis.h"
14#include "TF1.h"
15#include "TMath.h"
16
17
18/** \class TGraphQQ
19\ingroup BasicGraphics
20
21This class allows to draw quantile-quantile plots
22
23Plots can be drawn for 2 datasets or for a dataset and a theoretical
24distribution function
25
26## 2 datasets:
27 Quantile-quantile plots are used to determine whether 2 samples come from
28 the same distribution.
29 A qq-plot draws the quantiles of one dataset against the quantile of the
30 the other. The quantiles of the dataset with fewer entries are on Y axis,
31 with more entries - on X axis.
32 A straight line, going through 0.25 and 0.75 quantiles is also plotted
33 for reference. It represents a robust linear fit, not sensitive to the
34 extremes of the datasets.
35 If the datasets come from the same distribution, points of the plot should
36 fall approximately on the 45 degrees line. If they have the same
37 distribution function, but location or scale different parameters,
38 they should still fall on the straight line, but not the 45 degrees one.
39 The greater their departure from the straight line, the more evidence there
40 is, that the datasets come from different distributions.
41 The advantage of qq-plot is that it not only shows that the underlying
42 distributions are different, but, unlike the analytical methods, it also
43 gives information on the nature of this difference: heavier tails,
44 different location/scale, different shape, etc.
45
46 Some examples of qqplots of 2 datasets:
47
48\image html graf_graphqq1.png
49
50## 1 dataset:
51 Quantile-quantile plots are used to determine if the dataset comes from the
52 specified theoretical distribution, such as normal.
53 A qq-plot draws quantiles of the dataset against quantiles of the specified
54 theoretical distribution.
55 (NOTE, that density, not CDF should be specified)
56 A straight line, going through 0.25 and 0.75 quantiles can also be plotted
57 for reference. It represents a robust linear fit, not sensitive to the
58 extremes of the dataset.
59 As in the 2 datasets case, departures from straight line indicate departures
60 from the specified distribution.
61
62 "The correlation coefficient associated with the linear fit to the data
63 in the probability plot (qq plot in our case) is a measure of the
64 goodness of the fit.
65 Estimates of the location and scale parameters of the distribution
66 are given by the intercept and slope. Probability plots can be generated
67 for several competing distributions to see which provides the best fit,
68 and the probability plot generating the highest correlation coefficient
69 is the best choice since it generates the straightest probability plot."
70
71 From "Engineering statistic handbook",
72
73 http://www.itl.nist.gov/div898/handbook/eda/section3/probplot.htm
74
75 Example of a qq-plot of a dataset from N(3, 2) distribution and
76 TMath::Gaus(0, 1) theoretical function. Fitting parameters
77 are estimates of the distribution mean and sigma.
78
79\image html graf_graphqq2.png
80
81References:
82
83http://www.itl.nist.gov/div898/handbook/eda/section3/qqplot.htm
84
85http://www.itl.nist.gov/div898/handbook/eda/section3/probplot.htm
86*/
87
88////////////////////////////////////////////////////////////////////////////////
89/// default constructor
90
94
95////////////////////////////////////////////////////////////////////////////////
96/// Creates a quantile-quantile plot of dataset x.
97/// Theoretical distribution function can be defined later by SetFunction method
98
100 : TGraph(n)
101{
102 Int_t *index = new Int_t[n];
104 for (Int_t i=0; i<fNpoints; i++)
105 fY[i] = x[index[i]];
106 delete [] index;
107}
108
109////////////////////////////////////////////////////////////////////////////////
110/// Creates a quantile-quantile plot of dataset x against function f
111
113 : TGraph(n)
114{
115 fNy0 = 0;
116
117 Int_t *index = new Int_t[n];
119 for (Int_t i=0; i<fNpoints; i++)
120 fY[i] = x[index[i]];
121 delete [] index;
122 fF = f;
124}
125
126////////////////////////////////////////////////////////////////////////////////
127/// Creates a quantile-quantile plot of dataset x against dataset y
128/// Parameters nx and ny are respective array sizes
129
131{
132 fNpoints = (nx <= ny) ? nx : ny;
133
134 if (!CtorAllocate()) return;
135
136 Int_t *index = new Int_t[TMath::Max(nx, ny)];
138 if (nx <=ny){
139 for (Int_t i=0; i<fNpoints; i++)
140 fY[i] = x[index[i]];
142 if (nx==ny){
143 for (Int_t i=0; i<fNpoints; i++)
144 fX[i] = y[index[i]];
145 fY0 = nullptr;
146 Quartiles();
147 } else {
148 fNy0 = ny;
149 fY0 = new Double_t[ny];
150 for (Int_t i=0; i<ny; i++)
151 fY0[i] = y[i];
153 }
154 } else {
155 fNy0 = nx;
156 fY0 = new Double_t[nx];
157 for (Int_t i=0; i<nx; i++)
158 fY0[i] = x[index[i]];
160 for (Int_t i=0; i<ny; i++)
161 fY[i] = y[index[i]];
163 }
164
165 delete [] index;
166}
167
168////////////////////////////////////////////////////////////////////////////////
169/// Destroys a TGraphQQ
170
172{
173 if (fY0)
174 delete [] fY0;
175 if (fF)
176 fF = nullptr;
177}
178
179////////////////////////////////////////////////////////////////////////////////
180/// Computes quantiles of theoretical distribution function
181
183{
184 if (!fF) return;
185 TString s = fF->GetTitle();
186 Double_t pk;
187 if (s.Contains("TMath::Gaus") || s.Contains("gaus")){
188 //use plotting positions optimal for normal distribution
189 for (Int_t k=1; k<=fNpoints; k++){
190 pk = (k-0.375)/(fNpoints+0.25);
192 }
193 } else {
195 if (fNpoints > 10){
196 for (Int_t k=1; k<=fNpoints; k++)
197 prob[k-1] = (k-0.5)/fNpoints;
198 } else {
199 for (Int_t k=1; k<=fNpoints; k++)
200 prob[k-1] = (k-0.375)/(fNpoints+0.25);
201 }
202 //fF->GetQuantiles(fNpoints, prob, fX);
204 delete [] prob;
205 }
206
207 Quartiles();
208}
209
210////////////////////////////////////////////////////////////////////////////////
211/// When sample sizes are not equal, computes quantiles of the bigger sample
212/// by linear interpolation
213
215{
216 if (!fY0) return;
217
218 Double_t pi, pfrac;
219 Int_t pint;
220 for (Int_t i=0; i<fNpoints-1; i++){
221 pi = (fNy0-1)*Double_t(i)/Double_t(fNpoints-1);
223 pfrac = pi - pint;
224 fX[i] = (1-pfrac)*fY0[pint]+pfrac*fY0[pint+1];
225 }
226 fX[fNpoints-1]=fY0[fNy0-1];
227
228 Quartiles();
229}
230
231////////////////////////////////////////////////////////////////////////////////
232/// compute quartiles
233/// a quartile is a 25 per cent or 75 per cent quantile
234
236{
237 Double_t prob[]={0.25, 0.75};
238 Double_t x[2];
239 Double_t y[2];
241 if (fY0)
243 else if (fF) {
244 TString s = fF->GetTitle();
245 if (s.Contains("TMath::Gaus") || s.Contains("gaus")){
246 x[0] = TMath::NormQuantile(0.25);
247 x[1] = TMath::NormQuantile(0.75);
248 } else
249 fF->GetQuantiles(2, x, prob);
250 }
251 else
253
254 fXq1=x[0]; fXq2=x[1]; fYq1=y[0]; fYq2=y[1];
255}
256
257////////////////////////////////////////////////////////////////////////////////
258///Sets the theoretical distribution function (density!)
259///and computes its quantiles
260
262{
263 fF = f;
265}
#define f(i)
Definition RSha256.hxx:104
constexpr Bool_t kFALSE
Definition RtypesCore.h:108
double Double_t
Double 8 bytes.
Definition RtypesCore.h:73
constexpr Bool_t kTRUE
Definition RtypesCore.h:107
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t index
1-Dim function class
Definition TF1.h:182
virtual Int_t GetQuantiles(Int_t n, Double_t *xp, const Double_t *p)
Compute Quantiles for density distribution of this function.
Definition TF1.cxx:2019
void Quartiles()
compute quartiles a quartile is a 25 per cent or 75 per cent quantile
Definition TGraphQQ.cxx:235
TGraphQQ()
default constructor
Definition TGraphQQ.cxx:91
TF1 * fF
theoretical density function, if specified
Definition TGraphQQ.h:26
Double_t fYq1
y1 coordinate of the interquartile line
Definition TGraphQQ.h:23
Double_t * fY0
! second dataset, if specified
Definition TGraphQQ.h:25
void SetFunction(TF1 *f)
Sets the theoretical distribution function (density!) and computes its quantiles.
Definition TGraphQQ.cxx:261
~TGraphQQ() override
Destroys a TGraphQQ.
Definition TGraphQQ.cxx:171
Double_t fXq2
x2 coordinate of the interquartile line
Definition TGraphQQ.h:22
Int_t fNy0
size of the fY0 dataset
Definition TGraphQQ.h:20
void MakeFunctionQuantiles()
Computes quantiles of theoretical distribution function.
Definition TGraphQQ.cxx:182
Double_t fXq1
x1 coordinate of the interquartile line
Definition TGraphQQ.h:21
Double_t fYq2
y2 coordinate of the interquartile line
Definition TGraphQQ.h:24
void MakeQuantiles()
When sample sizes are not equal, computes quantiles of the bigger sample by linear interpolation.
Definition TGraphQQ.cxx:214
A TGraph is an object made of two arrays X and Y with npoints each.
Definition TGraph.h:41
Int_t fNpoints
Number of points <= fMaxSize.
Definition TGraph.h:46
Double_t * fY
[fNpoints] array of Y points
Definition TGraph.h:48
Bool_t CtorAllocate()
In constructors set fNpoints than call this method.
Definition TGraph.cxx:806
Double_t * fX
[fNpoints] array of X points
Definition TGraph.h:47
const char * GetTitle() const override
Returns title of object.
Definition TNamed.h:50
Basic string class.
Definition TString.h:138
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:640
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
const Int_t n
Definition legend1.C:16
Short_t Max(Short_t a, Short_t b)
Returns the largest of a and b.
Definition TMathBase.h:251
Int_t FloorNint(Double_t x)
Returns the nearest integer of TMath::Floor(x).
Definition TMath.h:697
Double_t NormQuantile(Double_t p)
Computes quantiles for standard normal distribution N(0, 1) at probability p.
Definition TMath.cxx:2459
void Quantiles(Int_t n, Int_t nprob, Double_t *x, Double_t *quantiles, Double_t *prob, Bool_t isSorted=kTRUE, Int_t *index=nullptr, Int_t type=7)
Computes sample quantiles, corresponding to the given probabilities.
Definition TMath.cxx:1207
void Sort(Index n, const Element *a, Index *index, Bool_t down=kTRUE)
Sort the n elements of the array a of generic templated type Element.
Definition TMathBase.h:432