0
$\begingroup$

So, i have a dataset which behaviour could be approached by the power law:

$$ N(x)= a - bx^{-1/c} $$

So i am doing this to get the regression: 

nlm = NonlinearModelFit[data, a - b*x^(-1/c), {{a, 0}, {b, 0}, {c,2}}, x, MaxIterations -> \[Infinity]]

So i get the following result:

enter image description here

Then i get these standard deviations for a, b, c:

enter image description here

Now i am using the gaussian law for error propagation to calculate the standard deviation of N(x):

\[Sqrt]((D[a - b*x^(-1/c), a]*11426.108699236138`)^2 + (D[a - b*x^(-1/c), b]*1.434690881054167`*^8)^2 + (D[a - b*x^(-1/c),c]*0.03626373539623972`)^2)

a = 3.2763042912759576`*^6;
b = 1.4018942844796562`*^9;
c = 2.320026692762461`;

Then i am getting the resulting function:

enter image description here

When i plot the resulting function sigma N(x), i am getting this graph:

enter image description here

So the error term is declines for x going to infinity.

Now my goal is to get the standard deviation for:

$$ x \rightarrow \infty $$

So i have done this:

enter image description here

So i am not that experienced in error calculations, and i am not sure if my action was right. I would be thankful if somebody could have a look on it.

|improve this question
$\endgroup$
  • $\begingroup$ With NonlinearModelFit you are fitting $N(x)=a−b x^{−1/c} + \epsilon$ where $\epsilon \sim N(0,\sigma^2)$. In other words, how you include the residual error is critical. Do the observed residuals support such an error structure? Are you looking for an estimate of the limit as $x\rightarrow \infty$ of $SE(\hat{a}-\hat{b} x^{-1/\hat{c}})$ (standard error of the mean prediction) or $SE(\hat{a}-\hat{b} x^{-1/\hat{c}}+\epsilon)$ (standard error of a single prediction)? $\endgroup$ – JimB Jan 27 '18 at 16:53
  • 1
    $\begingroup$ Once you have determined the appropriate error structure, you can have Mathematica do all of the work for you as it will calculate MeanPredictionBands and SinglePredictionBands. From those equations you can extract the associated standard errors. $\endgroup$ – JimB Jan 27 '18 at 17:10
  • $\begingroup$ @JimB Thank you, this helped! I think I am looking for for an estimate of the standard error of a single prediction. $\endgroup$ – Chopin Jan 28 '18 at 11:50
2
$\begingroup$

After thinking about this I think it's much simpler than I indicated in the comments.

The prediction for a single value of $x$ is $\hat{a}+\hat{b} x^{-1/\hat{c}}+\epsilon$. The limit of that prediction as $x\rightarrow \infty$ is $\hat{a}+\epsilon$ which you can check with the following command:

Limit of prediction

So you'd want to find the variance of $\hat{a}+\epsilon$ which is estimated to be

nlm["CovarianceMatrix"][[1, 1]] + nlm["EstimatedVariance"]

The more complicated approach I gave in the comments should get you the same number.

|improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.