2
$\begingroup$

Question:

I have some data, with a y=Ae^(-kt) model applied (It's a radioactive decay), thanks to answers on another question, the model is applied and returns the relevant numbers, but I need the value for the Chi^2 as a quantification for goodness of fit. Can anyone help me with how to about this? I've used PearsonChiSquaredFit in the example below.

Minimalised Example:

dataHist5 = {{Around[16.5, 1.5], Around[77.8, 8.8]}, {Around[34.5, 1.5], Around[60.5, 8.0]}, {Around[52.5, 1.5], Around[63.8, 8.0]}, {Around[106.5, 1.5], Around[42.4, 6.5]}, {Around[124.5, 1.5], Around[41.7, 6.5]}, {Around[142.5, 1.5], Around[14.6, 3.8]}, {Around[160.5, 1.5], Around[33.9, 5.8]}, {Around[178.5, 1.5], Around[29.4, 5.4]}}} data = Transpose[{dataHist5[[All, 1, 1]], dataHist5[[All, 2, 1]]}]; ListLogPlot[data] logData = data; logData[[All, 2]] = Log[data[[All, 2]]]/data[[All, 1]]; nlm = NonlinearModelFit[logData, loga/t - k, {loga, k}, t]; nlm["BestFitParameters"] Show[ListPlot[data], Plot[Exp[t nlm[t]], {t, 1, 700}]] ListPlot[Transpose[{logData[[All, 1]], nlm["FitResiduals"]}], PlotRange -> All] logData = data; logData[[All, 2]] = Log[data[[All, 2]]]/data[[All, 1]]; nlm = NonlinearModelFit[logData, loga/t - Log[2]/halfLife, {loga, halfLife}, t]; nlm["ParameterTable"] PearsonChiSquareTest[logData, nlm]

This returns:

PearsonChiSquareTest::rctnlndst: The argument FittedModel[-0.00676871+4.44614/t] at position 2 should be a valid distribution or a rectangular array of real numbers with length greater than the dimension of the array. The dimensionality of the arguments at positions 1 and 2 must match. {{16.5, 0.263887}, {34.5, 0.118917}, {52.5, 0.0791572}, {106.5, 0.0351845}, {124.5, 0.0299639}, {142.5, 0.0188142}, {160.5, 0.0219527}, {178.5, 0.0189411}, {196.5, 0.0178705}, {214.5, 0.0159942}, {232.5, 0.0147837}, {250.5, 0.0122477}, {268.5, 0.00543246}, {286.5, 0.00647923}, {322.5, 0.00624776}, {340.5, 0.00441726}, {358.5, 0.00668869}, {376.5, 0.00700945}, {394.5, 0.00668963}, {466.5, -0.00109502}, {484.5, 0.00162736}, {520.5, 0.00430492}, {538.5, 0.00262022}, {646.5, 0.00121958}, {682.5, -0.000748462}}, \!\(\* TagBox[ RowBox[{"FittedModel", "[", TagBox[ PanelBox[ TagBox[ RowBox[{ RowBox[{"-", "0.00676871191449669`"}], "+", FractionBox["4.4461438703116665`", "t"]}], Short[#, 2]& ], FrameMargins->5], Editable -> False], "]"}], InterpretTemplate[ FittedModel[{ "Nonlinear", {$CellContext`loga -> 4.4461438703116665`, \ $CellContext`halfLife -> 102.40459179174366`}, ...

Exponential Decay Graph of the expanded data set: Expanded data set e^-x graph

EDIT - R^2: I managed to find this for finding R^2. Ideally I would like to be able to return something like this for Chi^2

Grid[Transpose[{#, nlm[#]} &[{"AdjustedRSquared", "AIC", "BIC", "RSquared"}]], Alignment -> Left]

For me this returns: R^2 = 0.998555 Though that does look a bit high for the fit above.

EDIT - Real Data: Posting real data as requested:

dataHist5 = {{Around[16.5, 1.5], Around[77.8, 8.8]}, {Around[34.5, 1.5], Around[60.5, 8.0]}, {Around[52.5, 1.5], Around[63.8, 8.0]}, {Around[106.5, 1.5], Around[42.4, 6.5]}, {Around[124.5, 1.5], Around[41.7, 6.5]}, {Around[142.5, 1.5], Around[14.6, 3.8]}, {Around[160.5, 1.5], Around[33.9, 5.8]}, {Around[178.5, 1.5], Around[29.4, 5.4]}, {Around[196.5, 1.5], Around[33.5, 5.8]}, {Around[214.5, 1.5], Around[30.9, 5.6]}, {Around[232.5, 1.5], Around[31.1, 5.8]}, {Around[250.5, 1.5], Around[21.5, 4.6]}, {Around[268.5, 1.5], Around[4.3, 2.1]}, {Around[286.5, 1.5], Around[6.4, 2.5]}, {Around[322.5, 1.5], Around[7.5, 2.7]}, {Around[340.5, 1.5], Around[4.5, 2.1]}, {Around[358.5, 1.5], Around[11., 3.3]}, {Around[376.5, 1.5], Around[14.0, 3.7]}, {Around[394.5, 1.5], Around[14.0, 3.7]}, {Around[466.5, 1.5], Around[0.6, 0.7]}, {Around[502.5, 1.5], Around[2.2, 1.5]}, {Around[520.5, 1.5], Around[9.4, 3.1]}, {Around[538.5, 1.5], Around[4.1, 2.0]}, {Around[646.5, 1.5], Around[2.2, 1.5]}, {Around[682.5, 1.5], Around[0.6, 0.7]}}

EDIT - Recreating Answer

Show[ListPlot[dataHist5], Plot[fit[x], {x, 0, 800}, PlotRange -> All]] fit["ANOVATableSumsOfSquares"][[2]] (*chi^2*) fit["ANOVATableMeanSquares"][[2]] (*chi^2/dof*) uncertainties = dataHist5[[All, 2, 2]]; fit = NonlinearModelFit[rawDataHist5, A*Exp[-k*t], {A, k}, t, Weights -> 1/uncertainties^2];

Produces: Chi^2 = 90.68, Chi^2/DOF = 3.94. Which slightly differs from the quoted results of 90.75, 3.95 respectively. I'm not too sure why. The data isn't great so 3.95 seems like a realistic figure.

Improve this question
$\endgroup$
7
  • $\begingroup$ Is there some reason for needing a "Chi^2" statistic? The standard summary statistic (in most applied fields) is the root mean square error (RMSE) for assessing predictions and standard errors (or confidence intervals) for the regression parameters. $\endgroup$ Commented Mar 3, 2021 at 16:59
  • $\begingroup$ It was just suggested to me as a measure for goodness of fit, and matches previous times I've used it. But I am open to suggestions, how would I go about getting the RMSE in MMa? The ultimate goal is to get a figure for the half life, it's uncertainty, and a measure for how good the fit from the model was. It was also suggested from cross-validated that I try a Poisson GLM with log link model as opposed to standard nonlinearmodel, though I'm out my depth on statistics to know much there $\endgroup$ Commented Mar 3, 2021 at 17:10
  • $\begingroup$ The suggestion that you use a Poisson GLM was a good one. But those on CrossValidated assumed from your description that you had counts. And you need counts (i.e., integers) to use GeneralizedLinearModelFit. As a statistician when I want to do brain surgery, I consult a brain surgeon. Just a suggestion. $\endgroup$ Commented Mar 3, 2021 at 17:17
  • $\begingroup$ I think I might be missing the metaphor - is Cross Validated a poor place to obtain statistics advice, I was under the impression it was basically the Stats SE? $\endgroup$ Commented Mar 5, 2021 at 0:48
  • $\begingroup$ Sorry, my sarcasm was apparently very unclear. CrossValidated is great. What I was trying say was the you really ought to consult a statistician when you need statistical help. (And apparently not knowing that integer counts are required for Poisson regression seems to me a strong indication that you should really do so.) $\endgroup$ Commented Mar 5, 2021 at 0:51

1 Answer 1

Reset to default
3
$\begingroup$ 
nlm["ANOVATable"]

Will give you a summary of some variance metrics, including the Chi^2.

Under the Error row, you'll find

  • DF: degrees of freedom: 6
  • SS: Chi^2: 0.0000452471...
  • MS: Chi^2/dof: 7.54118*10^-6

If you want to extract the value programmatically, use

nlm["ANOVATableSumsOfSquares"][[2]]

PS:This is documented in the "Details and Options" tab of the NonlinearModelFitdocumentation.

For your real data edit:

see the following code:

rawDataHist5=dataHist5/.{Around[v_,__]:>v}; (*strip away the around's for the fit*) fit=NonlinearModelFit[rawDataHist5,A*Exp[-k*t],{A,k},t]; Show[ListPlot[dataHist5],Plot[fit[x],{x,0,800},PlotRange->All]] fit["ANOVATableSumsOfSquares"][[2]] (*chi^2*) fit["ANOVATableMeanSquares"][[2]] (*chi^2/dof*)

Which yields the plot Fit model with data

and the χ^2 of 1157.97 and a χ^2/dof of 50.35

However, your data can be weihted according to your uncertainties:

uncertainties=dataHist5[[All,2,2]]; fit=NonlinearModelFit[rawDataHist5,A*Exp[-k*t],{A,k},t,Weights->1/uncertainties^2];

In this case, the weights will be taken into account by calculating the χ^2 and we get χ^2 = 90.75 χ^2/dof = 3.95

Improve this answer
$\endgroup$
14
  • $\begingroup$ Thank you. Why was the pearsonchisquared command not working by the way? I'm getting the result of 0.000125, which seems like a very low figure, suggesting a big overfitting, but it doesn't look that way. And when I've fitted this graph on other things such as python and origin it's returned figures in the 10^2 range (which is still not great), do you know why there might be such a large difference? $\endgroup$ Commented Mar 2, 2021 at 16:12
  • $\begingroup$ If you look into the documentation for PearsonChiSquareTest, you see that the function does something totally different. Namely, using the Pearson X^2 test of checking for an underlying distribution of given data. Not evaluate the X^2 value for a fit. The X^2 value is only dependent on the fit result, so if origin and python came to the same parameter set, then all values should be the same. Maybe you compared your 'logData' fit with afit of the real data? If you post the real data this would be easy to show. $\endgroup$ Commented Mar 2, 2021 at 19:55
  • $\begingroup$ In the case of origin and python it's definitely not log data that's fitted because I haven't written any code that does that yet. I'd like to get the chi^2 of the exponential decay curve on the graph shown, and it's probanly safe to say I'm out my depth mathematica wise. I've included the full real data set as requested. $\endgroup$ Commented Mar 3, 2021 at 9:10
  • $\begingroup$ I edited the answer accordingly. The weighted answer of 3.95 for the chi^2/dof is realistically looking at your data. $\endgroup$ Commented Mar 3, 2021 at 10:21
  • $\begingroup$ Yeh, it's not great data is it? Might I ask why the weighting 1/uncertainty^2 is used rather than 1/uncertainty? And whether chi^2 or chi^2/dof is the most relevant metric here? $\endgroup$ Commented Mar 3, 2021 at 12:35

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.