Skip to content
Snippets Groups Projects
Verified Commit 29aa458b authored by Frank Sauerburger's avatar Frank Sauerburger
Browse files

Correct -1 in chi2 term

parent 2df8d2c9
No related branches found
No related tags found
No related merge requests found
Pipeline #1634 passed
......@@ -566,20 +566,20 @@ the uncertainties are the square roots of the expected number of events. If this
case, we have to compute the $`\chi^2`$ manually. The following example shows
both, the usage of `chisquare` and the manual computation. The print statements
for both methods produce the same output. Please note that
our model has five degrees of freedom since we have five free parameters. The
our model has five free parameters. The
total number of degrees of freedom is the number of data points reduced by the
number of degrees of freedom of our model.
of parameters of our model.
<!-- append decay.py -->
```python
# Degrees of freedom of our model.
dof = 5
# Delta degrees of freedom
ddof = 5
print("chi^2 from scipy:")
# Call chisquare and print.
chi2, p = scipy.stats.chisquare(count, fit_count, ddof=dof)
print(" chi2 / dof = %g / %d" % (chi2, len(count) - dof))
chi2, p = scipy.stats.chisquare(count, fit_count, ddof=ddof - 1)
print(" chi2 / dof = %g / %d" % (chi2, len(count) - ddof))
print(" p-value = %g" % p)
print() # print blank line
......@@ -591,17 +591,17 @@ uncertainty = np.sqrt(fit_count)
chi2 = ((count - fit_count)**2 / uncertainty**2).sum()
# Calculate p-value and print
p = scipy.stats.distributions.chi2.sf(chi2, len(count) - 1 - dof)
print(" chi2 / dof = %g / %d" % (chi2, len(count) - dof))
p = scipy.stats.distributions.chi2.sf(chi2, len(count) - ddof)
print(" chi2 / dof = %g / %d" % (chi2, len(count) - ddof))
print(" p-value = %g" % p)
```
The additional `-1` stems from the fact that we have a Poisson distribution
for which the uncertainties are derived from the number of events. The
test is not sensitive to the normalization, and therefore, we loose one degree
of freedom, cf. [Data Analysis in High Energy
Physics](https://www.wiley.com/en-us/Data+Analysis+in+High+Energy+Physics%3A+A+Practical+Guide+to+Statistical+Methods-p-9783527410583) by Behnke et at.
The additional `- 1` in the call of `chisquare` is necessary because the method
implicitly assumes that the total number of events is fixed. The method
therefore reduces the number of degrees of freedom by one. However, in our case,
the number of decays is not fixed and we have an additional degree
of freedom compared to what `chisquare` assumes. The `- 1` in the `ddof`
corrects for this.
If you run the `decay.py` you should see the following fit results.
<!-- console_output -->
......@@ -616,11 +616,11 @@ Optimal parameters:
chi^2 from scipy:
chi2 / dof = 120.577 / 123
p-value = 0.519418
p-value = 0.544943
Manual chi^2 test:
chi2 / dof = 120.577 / 123
p-value = 0.519418
p-value = 0.544943
```
Congratulations! You have mastered the first steps to analyze experimental data
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment