Logo ROOT  

Reference Guide
TestBinomial.C File Reference
Tutorials » Fit Tutorials

Detailed Description

View in nbviewer Open in SWAN Perform a fit to a set of data with binomial errors like those derived from the division of two histograms.

Three different fits are performed and compared:

  • simple least square fit to the divided histogram obtained from TH1::Divide with option b
  • least square fit to the TGraphAsymmErrors obtained from TGraphAsymmErrors::BayesDivide
  • likelihood fit performed on the dividing histograms using binomial statistics with the TBinomialEfficiency class

The first two methods are biased while the last one is statistical correct. Running the script passing an integer value n larger than 1, n fits are performed and the bias are also shown. To run the script :

to show the bias performing 100 fits for 1000 events per "experiment"

root[0]: .x TestBinomial.C+

to show the bias performing 100 fits for 1000 events per "experiment"

.x TestBinomial.C+(100, 1000)
32 68
****************************************
Minimizer is Minuit2
Chi2 = 0.350714
NDf = 3
Edm = 1.27001e-06
NCalls = 81
p0 = 0.694132 +/- 0.210029
p1 = 19.3471 +/- 5.85483
p2 = 5.2245 +/- 5.11013
****************************************
Minimizer is Minuit2
Chi2 = 1.64738
NDf = 5
Edm = 1.01171e-06
NCalls = 137
p0 = 0.677329 +/- 0.135567 (limited)
p1 = 15.5028 +/- 4.11118
p2 = 4.10441 +/- 2.77044
****************************************
Minimizer is Minuit2
Chi2 = 0.928456
NDf = 6
Edm = 1.46917e-08
NCalls = 75
p0 = 0.578806 +/- 0.128412 (limited)
p1 = 19.8914 +/- 2.83678
p2 = 3.39323 +/- 1.96156
****************************************
Minimizer is Minuit2
Chi2 = 2.56551
NDf = 7
Edm = 1.19556e-08
NCalls = 109
p0 = 0.57261 +/- 0.16552 (limited)
p1 = 21.4104 +/- 5.10611
p2 = 5.58786 +/- 3.30889
****************************************
Minimizer is Minuit2
Chi2 = 7.37095
NDf = 10
Edm = 1.54804e-06
NCalls = 87
p0 = 0.667689 +/- 0.152841 (limited)
p1 = 24.0146 +/- 2.40617
p2 = 3.63777 +/- 1.74672
****************************************
Minimizer is Minuit2
Chi2 = 0.747776
NDf = 6
Edm = 2.75121e-06
NCalls = 76
p0 = 0.651877 +/- 0.155245 (limited)
p1 = 18.9865 +/- 2.68321
p2 = 2.68989 +/- 1.28756
****************************************
Minimizer is Minuit2
Chi2 = 0.433039
NDf = 4
Edm = 3.07902e-07
NCalls = 104
p0 = 0.500967 +/- 0.254652 (limited)
p1 = 19.4199 +/- 9.30614
p2 = 6.46947 +/- 8.16624
****************************************
Minimizer is Minuit2
Chi2 = 4.49197
NDf = 9
Edm = 4.04441e-06
NCalls = 113
p0 = 0.623236 +/- 0.107197 (limited)
p1 = 20.8129 +/- 1.0507
p2 = 1.22179 +/- 1.09027
****************************************
Minimizer is Minuit2
Chi2 = 2.54966
NDf = 6
Edm = 1.48922e-07
NCalls = 93
p0 = 0.382687 +/- 0.226798 (limited)
p1 = 18.0539 +/- 6.38962
p2 = 3.0464 +/- 5.58751
****************************************
Minimizer is Minuit2
Chi2 = 1.46531
NDf = 7
Edm = 2.36369e-06
NCalls = 90
p0 = 0.575891 +/- 0.187613 (limited)
p1 = 17.6126 +/- 7.68638
p2 = 2.99185 +/- 6.18408
****************************************
Minimizer is Minuit2
Chi2 = 3.31728
NDf = 6
Edm = 8.18927e-07
NCalls = 64
p0 = 0.59853 +/- 0.124912 (limited)
p1 = 17.9048 +/- 3.49062
p2 = 3.52241 +/- 2.09942
****************************************
Minimizer is Minuit2
Chi2 = 0.413866
NDf = 6
Edm = 5.34616e-09
NCalls = 82
p0 = 0.654352 +/- 0.198145 (limited)
p1 = 20.6178 +/- 3.66715
p2 = 4.03019 +/- 2.20627
****************************************
Minimizer is Minuit2
Chi2 = 1.21002
NDf = 4
Edm = 2.90565e-08
NCalls = 82
p0 = 0.651685 +/- 0.159188 (limited)
p1 = 19.3369 +/- 2.66645
p2 = 2.81842 +/- 2.20413
****************************************
Minimizer is Minuit2
Chi2 = 0.752979
NDf = 7
Edm = 5.34806e-06
NCalls = 101
p0 = 0.643852 +/- 0.170146 (limited)
p1 = 21.01 +/- 4.32394
p2 = 5.36976 +/- 3.12567
****************************************
Minimizer is Minuit2
Chi2 = 2.5161
NDf = 8
Edm = 6.68524e-10
NCalls = 61
p0 = 0.628793 +/- 0.103749 (limited)
p1 = 17.4724 +/- 2.12101
p2 = 2.59351 +/- 1.42575
****************************************
Minimizer is Minuit2
Chi2 = 1.27043
NDf = 8
Edm = 6.23845e-06
NCalls = 72
p0 = 0.561127 +/- 0.138134 (limited)
p1 = 18.1601 +/- 3.22498
p2 = 3.50338 +/- 1.84312
****************************************
Minimizer is Minuit2
Chi2 = 2.84972
NDf = 8
Edm = 1.39524e-06
NCalls = 90
p0 = 0.693419 +/- 0.285373 (limited)
p1 = 24.7834 +/- 7.96871
p2 = 5.95505 +/- 3.88805
****************************************
Minimizer is Minuit2
Chi2 = 3.67016
NDf = 5
Edm = 1.80141e-06
NCalls = 82
p0 = 0.502645 +/- 0.176121 (limited)
p1 = 20.1539 +/- 2.79976
p2 = 2.70329 +/- 1.65981
****************************************
Invalid FitResult (status = 3 )
****************************************
Minimizer is Minuit2
Chi2 = 0.936088
NDf = 1
Edm = 0.106873
NCalls = 83
p0 = 0.41816 +/- 0.669502 (limited)
p1 = 18.3049 +/- 3.25491
p2 = 3.06774 +/- 40.715
****************************************
Minimizer is Minuit2
Chi2 = 3.77203
NDf = 8
Edm = 1.56279e-05
NCalls = 84
p0 = 0.642277 +/- 0.195607 (limited)
p1 = 22.9222 +/- 4.97062
p2 = 5.05681 +/- 2.69739
****************************************
Minimizer is Minuit2
Chi2 = 0.47406
NDf = 6
Edm = 1.7448e-06
NCalls = 88
p0 = 0.635811 +/- 0.173792 (limited)
p1 = 21.0478 +/- 3.90109
p2 = 3.91609 +/- 2.85302
****************************************
Minimizer is Minuit2
Chi2 = 2.39556
NDf = 8
Edm = 4.51695e-08
NCalls = 81
p0 = 0.564008 +/- 0.188334 (limited)
p1 = 19.2305 +/- 5.64705
p2 = 4.00431 +/- 3.2294
****************************************
Minimizer is Minuit2
Chi2 = 2.71864
NDf = 3
Edm = 2.02651e-07
NCalls = 84
p0 = 0.430955 +/- 0.103386 (limited)
p1 = 15.3556 +/- 1.58614
p2 = 1.18646 +/- 1.05279
****************************************
Minimizer is Minuit2
Chi2 = 3.91497
NDf = 3
Edm = 9.71484e-06
NCalls = 172
p0 = 0.999991 +/- 0.973989 (limited)
p1 = 26.3679 +/- 4.16357
p2 = 10.2085 +/- 7.97289
****************************************
Minimizer is Minuit2
Chi2 = 3.69474
NDf = 4
Edm = 1.0286e-07
NCalls = 124
p0 = 0.566565 +/- 0.300627 (limited)
p1 = 22.9958 +/- 7.83199
p2 = 6.59909 +/- 5.25595
****************************************
Minimizer is Minuit2
Chi2 = 1.01429
NDf = 10
Edm = 2.60252e-06
NCalls = 94
p0 = 0.610818 +/- 0.175926 (limited)
p1 = 23.3058 +/- 5.77903
p2 = 5.09666 +/- 2.59915
****************************************
Minimizer is Minuit2
Chi2 = 3.63097
NDf = 7
Edm = 2.77767e-07
NCalls = 123
p0 = 0.643028 +/- 0.225348 (limited)
p1 = 21.4306 +/- 4.46994
p2 = 3.95327 +/- 2.04398
****************************************
Minimizer is Minuit2
Chi2 = 1.7669
NDf = 7
Edm = 1.75691e-07
NCalls = 98
p0 = 0.73821 +/- 0.220641 (limited)
p1 = 23.671 +/- 5.18271
p2 = 5.53649 +/- 2.72352
****************************************
Minimizer is Minuit2
Chi2 = 2.70826
NDf = 8
Edm = 5.00596e-07
NCalls = 106
p0 = 0.662591 +/- 0.352545 (limited)
p1 = 24.2783 +/- 12.3384
p2 = 7.60846 +/- 6.66009
****************************************
Minimizer is Minuit2
Chi2 = 9.74855
NDf = 4
Edm = 3.46464e-07
NCalls = 124
p0 = 0.606488 +/- 0.231103 (limited)
p1 = 22.9147 +/- 5.11615
p2 = 4.88595 +/- 2.18611
****************************************
Minimizer is Minuit2
Chi2 = 3.10119
NDf = 5
Edm = 3.71723e-07
NCalls = 234
p0 = 0.999998 +/- 0.501156 (limited)
p1 = 28.5973 +/- 4.11085
p2 = 8.15931 +/- 3.01069
****************************************
Minimizer is Minuit2
Chi2 = 2.78024
NDf = 7
Edm = 1.5102e-06
NCalls = 82
p0 = 0.65308 +/- 0.206386 (limited)
p1 = 23.0571 +/- 4.95138
p2 = 4.86356 +/- 2.47873
****************************************
Minimizer is Minuit2
Chi2 = 3.05436
NDf = 5
Edm = 1.66356e-05
NCalls = 297
p0 = 0.999547 +/- 0.735933 (limited)
p1 = 47.3138 +/- 36.6961
p2 = 43.5419 +/- 61.6394
****************************************
Minimizer is Minuit2
Chi2 = 0.139804
NDf = 2
Edm = 1.06491e-06
NCalls = 91
p0 = 0.660157 +/- 0.258707 (limited)
p1 = 19.2592 +/- 7.6936
p2 = 6.44669 +/- 4.78428
****************************************
Minimizer is Minuit2
Chi2 = 2.55321
NDf = 7
Edm = 5.89829e-06
NCalls = 89
p0 = 0.643594 +/- 0.194962 (limited)
p1 = 21.9893 +/- 5.24123
p2 = 5.90243 +/- 2.90213
****************************************
Minimizer is Minuit2
Chi2 = 2.36128
NDf = 6
Edm = 6.99629e-07
NCalls = 88
p0 = 0.560279 +/- 0.13121 (limited)
p1 = 19.8526 +/- 3.03421
p2 = 2.47661 +/- 2.23047
****************************************
Minimizer is Minuit2
Chi2 = 1.33423
NDf = 7
Edm = 2.89896e-06
NCalls = 102
p0 = 0.564562 +/- 0.203287 (limited)
p1 = 25.8413 +/- 9.79483
p2 = 8.35761 +/- 6.07925
****************************************
Minimizer is Minuit2
Chi2 = 1.74047
NDf = 7
Edm = 7.04231e-07
NCalls = 135
p0 = 0.594764 +/- 0.187954 (limited)
p1 = 25.3877 +/- 8.57244
p2 = 6.52667 +/- 4.25157
****************************************
Minimizer is Minuit2
Chi2 = 0.310648
NDf = 7
Edm = 2.45348e-06
NCalls = 169
p0 = 0.539725 +/- 0.760718 (limited)
p1 = 12.9065 +/- 72.7676
p2 = 22.8828 +/- 100.306
****************************************
Minimizer is Minuit2
Chi2 = 2.52738
NDf = 10
Edm = 1.515e-08
NCalls = 119
p0 = 0.631088 +/- 0.21992 (limited)
p1 = 26.5773 +/- 9.03434
p2 = 8.9826 +/- 4.83446
****************************************
Minimizer is Minuit2
Chi2 = 2.51398
NDf = 6
Edm = 6.94471e-07
NCalls = 84
p0 = 0.484246 +/- 0.226988 (limited)
p1 = 18.2432 +/- 5.68591
p2 = 3.93586 +/- 3.23302
****************************************
Minimizer is Minuit2
Chi2 = 0.946468
NDf = 5
Edm = 5.45857e-06
NCalls = 74
p0 = 0.584615 +/- 0.133896 (limited)
p1 = 19.3728 +/- 3.81107
p2 = 4.20495 +/- 2.53748
****************************************
Minimizer is Minuit2
Chi2 = 5.46724
NDf = 8
Edm = 1.55331e-06
NCalls = 99
p0 = 0.763745 +/- 0.233374 (limited)
p1 = 26.3444 +/- 5.03855
p2 = 6.07972 +/- 2.87905
****************************************
Minimizer is Minuit2
Chi2 = 0.716283
NDf = 7
Edm = 1.89199e-06
NCalls = 85
p0 = 0.422768 +/- 0.179544 (limited)
p1 = 17.3878 +/- 6.43713
p2 = 3.05887 +/- 4.49587
****************************************
Minimizer is Minuit2
Chi2 = 6.40375
NDf = 8
Edm = 5.62321e-07
NCalls = 104
p0 = 0.507615 +/- 0.198193 (limited)
p1 = 22.2924 +/- 5.65071
p2 = 4.8117 +/- 2.67857
****************************************
Minimizer is Minuit2
Chi2 = 3.18171
NDf = 9
Edm = 5.5664e-07
NCalls = 83
p0 = 0.694498 +/- 0.134559 (limited)
p1 = 21.5399 +/- 2.86559
p2 = 3.46906 +/- 1.42287
****************************************
Minimizer is Minuit2
Chi2 = 3.13918
NDf = 6
Edm = 4.11187e-06
NCalls = 101
p0 = 0.725096 +/- 0.145408 (limited)
p1 = 25.2372 +/- 3.82188
p2 = 4.88315 +/- 1.80529
****************************************
Minimizer is Minuit2
Chi2 = 3.15774
NDf = 6
Edm = 2.32136e-06
NCalls = 109
p0 = 0.757722 +/- 0.770168 (limited)
p1 = 26.0106 +/- 13.1967
p2 = 7.01687 +/- 5.25991
****************************************
Minimizer is Minuit2
Chi2 = 1.90008
NDf = 6
Edm = 9.90083e-06
NCalls = 99
p0 = 0.534999 +/- 0.270262 (limited)
p1 = 24.0906 +/- 7.8809
p2 = 5.63935 +/- 3.24361
****************************************
Minimizer is Minuit2
Chi2 = 2.77855
NDf = 8
Edm = 8.13609e-08
NCalls = 85
p0 = 0.702887 +/- 0.199272 (limited)
p1 = 20.9159 +/- 3.68349
p2 = 4.41504 +/- 2.04434
****************************************
Minimizer is Minuit2
Chi2 = 4.37217
NDf = 9
Edm = 2.72083e-07
NCalls = 82
p0 = 0.620959 +/- 0.15843 (limited)
p1 = 21.0013 +/- 3.72087
p2 = 4.56895 +/- 2.37229
****************************************
Minimizer is Minuit2
Chi2 = 1.04047
NDf = 4
Edm = 5.64332e-06
NCalls = 102
p0 = 0.593644 +/- 0.174313 (limited)
p1 = 21.0676 +/- 4.04434
p2 = 3.35047 +/- 2.25854
****************************************
Minimizer is Minuit2
Chi2 = 3.95319
NDf = 10
Edm = 1.32825e-07
NCalls = 96
p0 = 0.574996 +/- 0.152016 (limited)
p1 = 20.2114 +/- 5.1893
p2 = 4.26461 +/- 2.57398
****************************************
Minimizer is Minuit2
Chi2 = 1.2121
NDf = 7
Edm = 3.20149e-06
NCalls = 91
p0 = 0.498388 +/- 0.194573 (limited)
p1 = 16.3792 +/- 10.7582
p2 = 4.04759 +/- 14.7395
****************************************
Minimizer is Minuit2
Chi2 = 1.15843
NDf = 4
Edm = 3.66535e-06
NCalls = 105
p0 = 0.570211 +/- 0.236898 (limited)
p1 = 17.2578 +/- 15.7293
p2 = 9.7494 +/- 16.3952
****************************************
Minimizer is Minuit2
Chi2 = 4.96865
NDf = 8
Edm = 1.88709e-08
NCalls = 80
p0 = 0.666513 +/- 0.140936 (limited)
p1 = 21.2569 +/- 2.26557
p2 = 3.48413 +/- 2.03653
****************************************
Minimizer is Minuit2
Chi2 = 2.88079
NDf = 4
Edm = 8.70733e-06
NCalls = 96
p0 = 0.803697 +/- 0.275939 (limited)
p1 = 22.1124 +/- 5.26388
p2 = 5.42701 +/- 3.06391
****************************************
Minimizer is Minuit2
Chi2 = 3.96392
NDf = 5
Edm = 7.79729e-08
NCalls = 100
p0 = 0.694889 +/- 0.235894 (limited)
p1 = 23.9686 +/- 5.40538
p2 = 5.07476 +/- 3.36335
****************************************
Minimizer is Minuit2
Chi2 = 0.852791
NDf = 5
Edm = 8.88345e-09
NCalls = 84
p0 = 0.310283 +/- 0.164695 (limited)
p1 = 15.0842 +/- 5.86864
p2 = 2.43104 +/- 3.51103
****************************************
Minimizer is Minuit2
Chi2 = 1.32394
NDf = 6
Edm = 6.16879e-07
NCalls = 86
p0 = 0.672106 +/- 0.201669 (limited)
p1 = 23.2133 +/- 5.10201
p2 = 5.04859 +/- 2.50327
****************************************
Minimizer is Minuit2
Chi2 = 2.48447
NDf = 9
Edm = 3.99523e-06
NCalls = 92
p0 = 0.469573 +/- 0.142912 (limited)
p1 = 18.7922 +/- 3.85884
p2 = 2.98044 +/- 2.31127
****************************************
Minimizer is Minuit2
Chi2 = 4.05841
NDf = 7
Edm = 3.35768e-06
NCalls = 93
p0 = 0.667164 +/- 0.148637 (limited)
p1 = 21.9757 +/- 2.55643
p2 = 3.34362 +/- 2.015
****************************************
Minimizer is Minuit2
Chi2 = 3.14114
NDf = 6
Edm = 1.74933e-08
NCalls = 114
p0 = 0.532234 +/- 0.174899 (limited)
p1 = 20.0178 +/- 3.71281
p2 = 3.05935 +/- 2.11129
****************************************
Minimizer is Minuit2
Chi2 = 2.88405
NDf = 5
Edm = 1.25252e-05
NCalls = 77
p0 = 0.698497 +/- 0.183511 (limited)
p1 = 21.6581 +/- 4.18422
p2 = 4.04207 +/- 1.56968
****************************************
Minimizer is Minuit2
Chi2 = 0.624725
NDf = 4
Edm = 6.43892e-06
NCalls = 71
p0 = 0.380579 +/- 0.134096 (limited)
p1 = 13.9323 +/- 2.56334
p2 = 1.96554 +/- 1.64705
****************************************
Minimizer is Minuit2
Chi2 = 1.67529
NDf = 7
Edm = 4.64318e-08
NCalls = 80
p0 = 0.616278 +/- 0.1576 (limited)
p1 = 18.812 +/- 2.91455
p2 = 3.21096 +/- 1.42859
****************************************
Minimizer is Minuit2
Chi2 = 0.374196
NDf = 2
Edm = 8.80185e-06
NCalls = 80
p0 = 0.685904 +/- 0.174442 (limited)
p1 = 19.998 +/- 1.62445
p2 = 1.71179 +/- 1.75684
****************************************
Minimizer is Minuit2
Chi2 = 2.46931
NDf = 7
Edm = 3.32581e-07
NCalls = 93
p0 = 0.640812 +/- 0.135875 (limited)
p1 = 22.3459 +/- 3.96589
p2 = 4.39764 +/- 2.51185
****************************************
Minimizer is Minuit2
Chi2 = 2.49631
NDf = 3
Edm = 6.6234e-06
NCalls = 102
p0 = 0.611704 +/- 0.35871 (limited)
p1 = 20.7057 +/- 8.42105
p2 = 5.08053 +/- 3.9281
****************************************
Minimizer is Minuit2
Chi2 = 0.418946
NDf = 2
Edm = 7.44661e-07
NCalls = 53
p0 = 0.60361 +/- 0.106733 (limited)
p1 = 15.7806 +/- 2
p2 = 0.105266 +/- 2
****************************************
Minimizer is Minuit2
Chi2 = 3.54456
NDf = 5
Edm = 2.33566e-06
NCalls = 85
p0 = 0.527397 +/- 0.217626 (limited)
p1 = 17.9928 +/- 6.38912
p2 = 5.4498 +/- 3.84602
****************************************
Minimizer is Minuit2
Chi2 = 3.23312
NDf = 8
Edm = 9.46056e-07
NCalls = 70
p0 = 0.515942 +/- 0.0980029 (limited)
p1 = 17.2963 +/- 1.44299
p2 = 1.9009 +/- 1.23849
****************************************
Minimizer is Minuit2
Chi2 = 5.28211
NDf = 11
Edm = 1.87776e-06
NCalls = 71
p0 = 0.612568 +/- 0.109907 (limited)
p1 = 19.3198 +/- 2.18618
p2 = 3.37801 +/- 1.31028
****************************************
Minimizer is Minuit2
Chi2 = 2.51155
NDf = 7
Edm = 2.60724e-06
NCalls = 74
p0 = 0.647126 +/- 0.147452 (limited)
p1 = 20.8474 +/- 2.85119
p2 = 3.27715 +/- 1.66359
****************************************
Minimizer is Minuit2
Chi2 = 2.33728
NDf = 6
Edm = 7.37813e-07
NCalls = 109
p0 = 0.624689 +/- 0.141987 (limited)
p1 = 19.3681 +/- 3.57572
p2 = 3.40927 +/- 2.17685
****************************************
Minimizer is Minuit2
Chi2 = 3.75642
NDf = 8
Edm = 5.68096e-06
NCalls = 120
p0 = 0.387014 +/- 0.0702296 (limited)
p1 = 11.5811 +/- 8.22684
p2 = 0.313997 +/- 4.44032
****************************************
Minimizer is Minuit2
Chi2 = 1.68054
NDf = 6
Edm = 5.00607e-06
NCalls = 79
p0 = 0.492763 +/- 0.200749 (limited)
p1 = 18.8758 +/- 5.209
p2 = 3.88829 +/- 2.85902
****************************************
Minimizer is Minuit2
Chi2 = 1.37515
NDf = 4
Edm = 3.40318e-06
NCalls = 106
p0 = 0.759679 +/- 0.234791 (limited)
p1 = 24.2145 +/- 5.44887
p2 = 4.93415 +/- 2.65778
****************************************
Minimizer is Minuit2
Chi2 = 2.35692
NDf = 4
Edm = 1.08773e-07
NCalls = 81
p0 = 0.373315 +/- 0.215845 (limited)
p1 = 12.4572 +/- 7.57085
p2 = 4.00427 +/- 8.97995
****************************************
Minimizer is Minuit2
Chi2 = 1.93137
NDf = 7
Edm = 4.79943e-09
NCalls = 113
p0 = 0.632872 +/- 0.231351 (limited)
p1 = 27.1811 +/- 8.5188
p2 = 6.83092 +/- 3.30519
****************************************
Minimizer is Minuit2
Chi2 = 3.28016
NDf = 6
Edm = 6.04893e-06
NCalls = 99
p0 = 0.696957 +/- 0.193092 (limited)
p1 = 22.9427 +/- 3.44124
p2 = 3.63739 +/- 2.36972
****************************************
Minimizer is Minuit2
Chi2 = 2.56902
NDf = 7
Edm = 5.20615e-07
NCalls = 94
p0 = 0.68832 +/- 0.19603 (limited)
p1 = 24.2988 +/- 4.41739
p2 = 4.93562 +/- 3.1133
****************************************
Minimizer is Minuit2
Chi2 = 3.04459
NDf = 7
Edm = 7.4353e-07
NCalls = 83
p0 = 0.659997 +/- 0.249426 (limited)
p1 = 22.9747 +/- 5.95617
p2 = 5.36022 +/- 2.66473
****************************************
Minimizer is Minuit2
Chi2 = 6.27397
NDf = 8
Edm = 5.29961e-07
NCalls = 81
p0 = 0.662433 +/- 0.15714 (limited)
p1 = 21.5909 +/- 3.29623
p2 = 4.09179 +/- 1.74745
****************************************
Minimizer is Minuit2
Chi2 = 1.94311
NDf = 9
Edm = 1.89133e-07
NCalls = 106
p0 = 0.575157 +/- 0.298322 (limited)
p1 = 25.6868 +/- 11.7079
p2 = 7.54114 +/- 5.74853
****************************************
Minimizer is Minuit2
Chi2 = 2.75509
NDf = 9
Edm = 9.98946e-08
NCalls = 69
p0 = 0.593401 +/- 0.101469 (limited)
p1 = 18.8032 +/- 2.57519
p2 = 2.53183 +/- 1.7654
****************************************
Minimizer is Minuit2
Chi2 = 2.54952
NDf = 6
Edm = 4.57669e-08
NCalls = 96
p0 = 0.660242 +/- 0.259374 (limited)
p1 = 25.0771 +/- 5.98352
p2 = 5.39543 +/- 2.98249
****************************************
Minimizer is Minuit2
Chi2 = 4.62283
NDf = 10
Edm = 1.36401e-08
NCalls = 89
p0 = 0.592713 +/- 0.183874 (limited)
p1 = 20.9349 +/- 5.53213
p2 = 6.39618 +/- 4.21336
****************************************
Minimizer is Minuit2
Chi2 = 5.3127
NDf = 6
Edm = 2.53153e-06
NCalls = 119
p0 = 0.534666 +/- 0.0847463 (limited)
p1 = 12.0856 +/- 13.9766
p2 = 0.685012 +/- 8.82276
****************************************
Minimizer is Minuit2
Chi2 = 0.968364
NDf = 7
Edm = 1.15106e-06
NCalls = 74
p0 = 0.669971 +/- 0.134042 (limited)
p1 = 19.9248 +/- 3.66311
p2 = 4.7186 +/- 3.31403
****************************************
Minimizer is Minuit2
Chi2 = 2.88345
NDf = 10
Edm = 2.24225e-06
NCalls = 102
p0 = 0.629842 +/- 0.174064 (limited)
p1 = 21.9766 +/- 6.0423
p2 = 6.83408 +/- 3.9968
****************************************
Minimizer is Minuit2
Chi2 = 1.80215
NDf = 6
Edm = 2.68316e-07
NCalls = 87
p0 = 0.495591 +/- 0.211145 (limited)
p1 = 19.645 +/- 5.2052
p2 = 3.59211 +/- 2.51889
****************************************
Minimizer is Minuit2
Chi2 = 3.79465
NDf = 6
Edm = 4.35871e-08
NCalls = 69
p0 = 0.611171 +/- 0.128277 (limited)
p1 = 17.467 +/- 2.93532
p2 = 3.02848 +/- 1.64432
****************************************
Minimizer is Minuit2
Chi2 = 2.748
NDf = 7
Edm = 5.78153e-09
NCalls = 140
p0 = 1 +/- 0.733966 (limited)
p1 = 28.7334 +/- 2.38732
p2 = 6.33261 +/- 1.85434
****************************************
Minimizer is Minuit2
Chi2 = 4.29305
NDf = 9
Edm = 3.06585e-06
NCalls = 82
p0 = 0.687721 +/- 0.16847 (limited)
p1 = 22.2779 +/- 3.60769
p2 = 4.13806 +/- 1.75728
****************************************
Minimizer is Minuit2
Chi2 = 1.74986
NDf = 2
Edm = 2.69957e-06
NCalls = 256
p0 = 0.999951 +/- 0.988093 (limited)
p1 = 33.6961 +/- 9.16268
p2 = 13.7955 +/- 8.24703
****************************************
Minimizer is Minuit2
Chi2 = 4.72905
NDf = 5
Edm = 2.89445e-07
NCalls = 67
p0 = 0.417412 +/- 0.107227 (limited)
p1 = 14.5796 +/- 2.95636
p2 = 2.22066 +/- 2.12966
****************************************
Minimizer is Minuit2
Chi2 = 2.45643
NDf = 6
Edm = 2.63752e-06
NCalls = 99
p0 = 0.554001 +/- 0.226655 (limited)
p1 = 23.9121 +/- 6.69531
p2 = 6.04267 +/- 3.44117
****************************************
Minimizer is Minuit2
Chi2 = 2.91637
NDf = 4
Edm = 1.377e-05
NCalls = 84
p0 = 0.769946 +/- 0.185015 (limited)
p1 = 20.5779 +/- 2.52006
p2 = 3.50084 +/- 1.35396
****************************************
Minimizer is Minuit2
Chi2 = 1.21412
NDf = 6
Edm = 1.16053e-08
NCalls = 141
p0 = 0.632519 +/- 0.106789 (limited)
p1 = 17.4834 +/- 1.66297
p2 = 1.44082 +/- 1.10945
#include "TBinomialEfficiencyFitter.h"
#include "TVirtualFitter.h"
#include "TH1.h"
#include "TRandom3.h"
#include "TF1.h"
#include "TFitResult.h"
#include "TStyle.h"
#include "TCanvas.h"
#include "TLegend.h"
#include "TPaveStats.h"
#include "Math/IntegratorOptions.h"
#include <cassert>
#include <iostream>
void TestBinomial(int nloop = 100, int nevts = 100, bool plot = false, bool debug = false, int seed = 111)
{
gStyle->SetMarkerStyle(20);
gStyle->SetLineWidth(2.0);
gStyle->SetOptStat(11);
TObjArray hbiasNorm;
hbiasNorm.Add(new TH1D("h0Norm", "Bias Histogram fit",100,-5,5));
hbiasNorm.Add(new TH1D("h1Norm","Bias Binomial fit",100,-5,5));
TObjArray hbiasThreshold;
hbiasThreshold.Add(new TH1D("h0Threshold", "Bias Histogram fit",100,-5,5));
hbiasThreshold.Add(new TH1D("h1Threshold","Bias Binomial fit",100,-5,5));
TObjArray hbiasWidth;
hbiasWidth.Add(new TH1D("h0Width", "Bias Histogram fit",100,-5,5));
hbiasWidth.Add(new TH1D("h1Width","Bias Binomial fit",100,-5,5));
TH1D* hChisquared = new TH1D("hChisquared",
"#chi^{2} probability (Baker-Cousins)", 200, 0.0, 1.0);
TVirtualFitter::SetDefaultFitter("Minuit2");
ROOT::Math::IntegratorOneDimOptions::SetDefaultIntegrator("Gauss");
// Note: in order to be able to use TH1::FillRandom() to generate
// pseudo-experiments, we use a trick: generate "selected"
// and "non-selected" samples independently. These are
// statistically independent and therefore can be safely
// added to yield the "before selection" sample.
// Define (arbitrarily?) a distribution of input events.
// Here: assume a x^(-2) distribution. Boundaries: [10, 100].
Double_t xmin =10, xmax = 100;
TH1D* hM2D = new TH1D("hM2D", "x^(-2) denominator distribution",
45, xmin, xmax);
TH1D* hM2N = new TH1D("hM2N", "x^(-2) numerator distribution",
45, xmin, xmax);
TH1D* hM2E = new TH1D("hM2E", "x^(-2) efficiency",
45, xmin, xmax);
TF1* fM2D = new TF1("fM2D", "(1-[0]/(1+exp(([1]-x)/[2])))/(x*x)",
xmin, xmax);
TF1* fM2N = new TF1("fM2N", "[0]/(1+exp(([1]-x)/[2]))/(x*x)",
xmin, xmax);
TF1* fM2Fit = new TF1("fM2Fit", "[0]/(1+exp(([1]-x)/[2]))",
xmin, xmax);
TF1* fM2Fit2 = 0;
TRandom3 rb(seed);
// First try: use a single set of parameters.
// For each try, we need to find the overall normalization
Double_t normalization = 0.80;
Double_t threshold = 25.0;
Double_t width = 5.0;
fM2D->SetParameter(0, normalization);
fM2D->SetParameter(1, threshold);
fM2D->SetParameter(2, width);
fM2N->SetParameter(0, normalization);
fM2N->SetParameter(1, threshold);
fM2N->SetParameter(2, width);
Double_t integralN = fM2N->Integral(xmin, xmax);
Double_t integralD = fM2D->Integral(xmin, xmax);
Double_t fracN = integralN/(integralN+integralD);
Int_t nevtsN = rb.Binomial(nevts, fracN);
Int_t nevtsD = nevts - nevtsN;
std::cout << nevtsN << " " << nevtsD << std::endl;
gStyle->SetOptFit(1111);
// generate many times to see the bias
for (int iloop = 0; iloop < nloop; ++iloop) {
// generate pseudo-experiments
hM2D->Reset();
hM2N->Reset();
hM2D->FillRandom(fM2D->GetName(), nevtsD);
hM2N->FillRandom(fM2N->GetName(), nevtsN);
hM2D->Add(hM2N);
// construct the "efficiency" histogram
hM2N->Sumw2();
hM2E->Divide(hM2N, hM2D, 1, 1, "b");
// Fit twice, using the same fit function.
// In the first (standard) fit, initialize to (arbitrary) values.
// In the second fit, use the results from the first fit (this
// makes it easier for the fit -- but the purpose here is not to
// show how easy or difficult it is to obtain results, but whether
// the CORRECT results are obtained or not!).
fM2Fit->SetParameter(0, 0.5);
fM2Fit->SetParameter(1, 15.0);
fM2Fit->SetParameter(2, 2.0);
fM2Fit->SetParError(0, 0.1);
fM2Fit->SetParError(1, 1.0);
fM2Fit->SetParError(2, 0.2);
TH1 * hf = fM2Fit->GetHistogram();
// std::cout << "Function values " << std::endl;
// for (int i = 1; i <= hf->GetNbinsX(); ++i)
// std::cout << hf->GetBinContent(i) << " ";
// std::cout << std::endl;
TCanvas* cEvt;
if (plot) {
cEvt = new TCanvas(Form("cEnv%d",iloop),
Form("plots for experiment %d", iloop),
1000, 600);
cEvt->Divide(1,2);
cEvt->cd(1);
hM2D->DrawCopy("HIST");
hM2N->SetLineColor(kRed);
hM2N->DrawCopy("HIST SAME");
cEvt->cd(2);
}
for (int fit = 0; fit < 2; ++fit) {
Int_t status = 0;
switch (fit) {
case 0:
{
// TVirtualPad * pad = gPad;
// new TCanvas();
// fM2Fit->Draw();
// gPad = pad;
TString optFit = "RN";
if (debug) optFit += TString("SV");
TFitResultPtr res = hM2E->Fit(fM2Fit, optFit);
if (plot) {
hM2E->DrawCopy("E");
fM2Fit->SetLineColor(kBlue);
fM2Fit->DrawCopy("SAME");
}
if (debug) res->Print();
status = res;
break;
}
case 1:
{
// if (fM2Fit2) delete fM2Fit2;
// fM2Fit2 = dynamic_cast<TF1*>(fM2Fit->Clone("fM2Fit2"));
fM2Fit2 = fM2Fit; // do not clone/copy the function
if (fM2Fit2->GetParameter(0) >= 1.0)
fM2Fit2->SetParameter(0, 0.95);
fM2Fit2->SetParLimits(0, 0.0, 1.0);
// TVirtualPad * pad = gPad;
// new TCanvas();
// fM2Fit2->Draw();
// gPad = pad;
TBinomialEfficiencyFitter bef(hM2N, hM2D);
TString optFit = "RI S";
if (debug) optFit += TString("V");
TFitResultPtr res = bef.Fit(fM2Fit2,optFit);
status = res;
if (status !=0) {
std::cerr << "Error performing binomial efficiency fit, result = "
<< status << std::endl;
res->Print();
continue;
}
if (plot) {
fM2Fit2->SetLineColor(kRed);
fM2Fit2->DrawCopy("SAME");
bool confint = (status == 0);
if (confint) {
// compute confidence interval on fitted function
auto htemp = fM2Fit2->GetHistogram();
ROOT::Fit::BinData xdata;
ROOT::Fit::FillData(xdata, fM2Fit2->GetHistogram() );
TGraphErrors gr(fM2Fit2->GetHistogram() );
res->GetConfidenceIntervals(xdata, gr.GetEY(), 0.68, false);
gr.SetFillColor(6);
gr.SetFillStyle(3005);
gr.DrawClone("4 same");
}
}
if (debug) {
res->Print();
}
}
}
if (status != 0) break;
Double_t fnorm = fM2Fit->GetParameter(0);
Double_t enorm = fM2Fit->GetParError(0);
Double_t fthreshold = fM2Fit->GetParameter(1);
Double_t ethreshold = fM2Fit->GetParError(1);
Double_t fwidth = fM2Fit->GetParameter(2);
Double_t ewidth = fM2Fit->GetParError(2);
if (fit == 1) {
fnorm = fM2Fit2->GetParameter(0);
enorm = fM2Fit2->GetParError(0);
fthreshold = fM2Fit2->GetParameter(1);
ethreshold = fM2Fit2->GetParError(1);
fwidth = fM2Fit2->GetParameter(2);
ewidth = fM2Fit2->GetParError(2);
hChisquared->Fill(fM2Fit2->GetProb());
}
TH1D* h = dynamic_cast<TH1D*>(hbiasNorm[fit]);
h->Fill((fnorm-normalization)/enorm);
h = dynamic_cast<TH1D*>(hbiasThreshold[fit]);
h->Fill((fthreshold-threshold)/ethreshold);
h = dynamic_cast<TH1D*>(hbiasWidth[fit]);
h->Fill((fwidth-width)/ewidth);
}
}
TCanvas* c1 = new TCanvas("c1",
"Efficiency fit biases",10,10,1000,800);
c1->Divide(2,2);
TH1D *h0, *h1;
c1->cd(1);
h0 = dynamic_cast<TH1D*>(hbiasNorm[0]);
h0->Draw("HIST");
h1 = dynamic_cast<TH1D*>(hbiasNorm[1]);
h1->SetLineColor(kRed);
h1->Draw("HIST SAMES");
TLegend* l1 = new TLegend(0.1, 0.75, 0.5, 0.9,
"plateau parameter", "ndc");
l1->AddEntry(h0, Form("histogram: mean = %4.2f RMS = \
%4.2f", h0->GetMean(), h0->GetRMS()), "l");
l1->AddEntry(h1, Form("binomial : mean = %4.2f RMS = \
%4.2f", h1->GetMean(), h1->GetRMS()), "l");
l1->Draw();
c1->cd(2);
h0 = dynamic_cast<TH1D*>(hbiasThreshold[0]);
h0->Draw("HIST");
h1 = dynamic_cast<TH1D*>(hbiasThreshold[1]);
h1->SetLineColor(kRed);
h1->Draw("HIST SAMES");
TLegend* l2 = new TLegend(0.1, 0.75, 0.5, 0.9,
"threshold parameter", "ndc");
l2->AddEntry(h0, Form("histogram: mean = %4.2f RMS = \
%4.2f", h0->GetMean(), h0->GetRMS()), "l");
l2->AddEntry(h1, Form("binomial : mean = %4.2f RMS = \
%4.2f", h1->GetMean(), h1->GetRMS()), "l");
l2->Draw();
c1->cd(3);
h0 = dynamic_cast<TH1D*>(hbiasWidth[0]);
h0->Draw("HIST");
h1 = dynamic_cast<TH1D*>(hbiasWidth[1]);
h1->SetLineColor(kRed);
h1->Draw("HIST SAMES");
TLegend* l3 = new TLegend(0.1, 0.75, 0.5, 0.9, "width parameter", "ndc");
l3->AddEntry(h0, Form("histogram: mean = %4.2f RMS = \
%4.2f", h0->GetMean(), h0->GetRMS()), "l");
l3->AddEntry(h1, Form("binomial : mean = %4.2f RMS = \
%4.2f", h1->GetMean(), h1->GetRMS()), "l");
l3->Draw();
c1->cd(4);
hChisquared->Draw("HIST");
}
int main() {
TestBinomial();
}
h
#define h(i)
Definition: RSha256.hxx:106
Int_t
int Int_t
Definition: RtypesCore.h:41
Double_t
double Double_t
Definition: RtypesCore.h:55
kRed
@ kRed
Definition: Rtypes.h:64
kBlue
@ kBlue
Definition: Rtypes.h:64
width
include TDocParser_001 C image html pict1_TDocParser_001 png width
Definition: TDocParser.cxx:121
xmin
float xmin
Definition: THbookFile.cxx:93
xmax
float xmax
Definition: THbookFile.cxx:93
Form
char * Form(const char *fmt,...)
gStyle
R__EXTERN TStyle * gStyle
Definition: TStyle.h:407
ROOT::Fit::BinData
Class describing the binned data sets : vectors of x coordinates, y values and optionally error on y ...
Definition: BinData.h:53
ROOT::Fit::FitResult::GetConfidenceIntervals
void GetConfidenceIntervals(unsigned int n, unsigned int stride1, unsigned int stride2, const double *x, double *ci, double cl=0.95, bool norm=false) const
get confidence intervals for an array of n points x.
Definition: FitResult.cxx:553
ROOT::Math::IntegratorOneDimOptions::SetDefaultIntegrator
static void SetDefaultIntegrator(const char *name)
Definition: IntegratorOptions.cxx:223
TAttFill::SetFillColor
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition: TAttFill.h:37
TAttFill::SetFillStyle
virtual void SetFillStyle(Style_t fstyle)
Set the fill area style.
Definition: TAttFill.h:39
TAttLine::SetLineWidth
virtual void SetLineWidth(Width_t lwidth)
Set the line width.
Definition: TAttLine.h:43
TAttLine::SetLineColor
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition: TAttLine.h:40
TAttMarker::SetMarkerStyle
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition: TAttMarker.h:40
TBinomialEfficiencyFitter
Binomial fitter for the division of two histograms.
Definition: TBinomialEfficiencyFitter.h:42
TCanvas
The Canvas class.
Definition: TCanvas.h:31
TCanvas::cd
TVirtualPad * cd(Int_t subpadnumber=0)
Set current canvas & pad.
Definition: TCanvas.cxx:696
TF1
1-Dim function class
Definition: TF1.h:211
TF1::SetParError
virtual void SetParError(Int_t ipar, Double_t error)
Set error for parameter number ipar.
Definition: TF1.cxx:3472
TF1::GetHistogram
virtual TH1 * GetHistogram() const
Return a pointer to the histogram used to visualise the function.
Definition: TF1.cxx:1564
TF1::GetParError
virtual Double_t GetParError(Int_t ipar) const
Return value of parameter number ipar.
Definition: TF1.cxx:1910
TF1::GetProb
virtual Double_t GetProb() const
Return the fit probability.
Definition: TF1.cxx:1935
TF1::Integral
virtual Double_t Integral(Double_t a, Double_t b, Double_t epsrel=1.e-12)
IntegralOneDim or analytical integral.
Definition: TF1.cxx:2502
TF1::SetParLimits
virtual void SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax)
Set limits for parameter ipar.
Definition: TF1.cxx:3497
TF1::DrawCopy
virtual TF1 * DrawCopy(Option_t *option="") const
Draw a copy of this function with its current attributes.
Definition: TF1.cxx:1347
TF1::SetParameter
virtual void SetParameter(Int_t param, Double_t value)
Definition: TF1.h:628
TF1::GetParameter
virtual Double_t GetParameter(Int_t ipar) const
Definition: TF1.h:506
TFitResultPtr
Provides an indirection to the TFitResult class and with a semantics identical to a TFitResult pointe...
Definition: TFitResultPtr.h:31
TFitResult::Print
virtual void Print(Option_t *option="") const
Print result of the fit, by default chi2, parameter values and errors.
Definition: TFitResult.cxx:44
TGraphErrors
A TGraphErrors is a TGraph with error bars.
Definition: TGraphErrors.h:26
TGraphErrors::GetEY
Double_t * GetEY() const
Definition: TGraphErrors.h:67
TH1D
1-D histogram with a double per channel (see TH1 documentation)}
Definition: TH1.h:614
TH1D::Reset
virtual void Reset(Option_t *option="")
Reset.
Definition: TH1.cxx:9729
TH1
The TH1 histogram class.
Definition: TH1.h:56
TH1::FillRandom
virtual void FillRandom(const char *fname, Int_t ntimes=5000)
Fill histogram following distribution in function fname.
Definition: TH1.cxx:3445
TH1::GetMean
virtual Double_t GetMean(Int_t axis=1) const
For axis = 1,2 or 3 returns the mean value of the histogram along X,Y or Z axis.
Definition: TH1.cxx:7074
TH1::Fit
virtual TFitResultPtr Fit(const char *formula, Option_t *option="", Option_t *goption="", Double_t xmin=0, Double_t xmax=0)
Fit histogram with function fname.
Definition: TH1.cxx:3808
TH1::Add
virtual Bool_t Add(TF1 *h1, Double_t c1=1, Option_t *option="")
Performs the operation: this = this + c1*f1 if errors are defined (see TH1::Sumw2),...
Definition: TH1.cxx:778
TH1::Fill
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3275
TH1::GetRMS
Double_t GetRMS(Int_t axis=1) const
Definition: TH1.h:311
TH1::DrawCopy
virtual TH1 * DrawCopy(Option_t *option="", const char *name_postfix="_copy") const
Copy this histogram and Draw in the current pad.
Definition: TH1.cxx:3045
TH1::Draw
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2998
TH1::Divide
virtual Bool_t Divide(TF1 *f1, Double_t c1=1)
Performs the operation: this = this/(c1*f1) if errors are defined (see TH1::Sumw2),...
Definition: TH1.cxx:2753
TH1::Sumw2
virtual void Sumw2(Bool_t flag=kTRUE)
Create structure to store sum of squares of weights.
Definition: TH1.cxx:8464
TLegend
This class displays a legend box (TPaveText) containing several legend entries.
Definition: TLegend.h:23
TLegend::AddEntry
TLegendEntry * AddEntry(const TObject *obj, const char *label="", Option_t *option="lpf")
Add a new entry to this legend.
Definition: TLegend.cxx:330
TLegend::Draw
virtual void Draw(Option_t *option="")
Draw this legend with its current attributes.
Definition: TLegend.cxx:423
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
TObjArray
An array of TObjects.
Definition: TObjArray.h:37
TObjArray::Add
void Add(TObject *obj)
Definition: TObjArray.h:74
TObject::DrawClone
virtual TObject * DrawClone(Option_t *option="") const
Draw a clone of this object in the current selected pad for instance with: gROOT->SetSelectedPad(gPad...
Definition: TObject.cxx:219
TPad::Divide
virtual void Divide(Int_t nx=1, Int_t ny=1, Float_t xmargin=0.01, Float_t ymargin=0.01, Int_t color=0)
Automatic pad generation by division.
Definition: TPad.cxx:1163
TRandom3
Random number generator class based on M.
Definition: TRandom3.h:27
TString
Basic string class.
Definition: TString.h:131
TStyle::SetOptStat
void SetOptStat(Int_t stat=1)
The type of information printed in the histogram statistics box can be selected via the parameter mod...
Definition: TStyle.cxx:1450
TStyle::SetOptFit
void SetOptFit(Int_t fit=1)
The type of information about fit parameters printed in the histogram statistics box can be selected ...
Definition: TStyle.cxx:1402
TVirtualFitter::SetDefaultFitter
static void SetDefaultFitter(const char *name="")
static: set name of default fitter
Definition: TVirtualFitter.cxx:245
main
int main(int argc, char **argv)
Definition: histspeedtest.cxx:752
c1
return c1
Definition: legend1.C:41
gr
TGraphErrors * gr
Definition: legend1.C:25
h1
TH1F * h1
Definition: legend1.C:5
ROOT::Fit::FillData
void FillData(BinData &dv, const TH1 *hist, TF1 *func=0)
fill the data vector from a TH1.
Definition: HFitInterface.cxx:107
Author
Rene Brun

Definition in file TestBinomial.C.