# Calculating standard deviation of power law regression going to infinity

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: Then i get these standard deviations for a, b, c: 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: When i plot the resulting function sigma N(x), i am getting this graph: 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: 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.

• 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)?
– JimB
Jan 27, 2018 at 16:53
• 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.
– JimB
Jan 27, 2018 at 17:10
• @JimB Thank you, this helped! I think I am looking for for an estimate of the standard error of a single prediction. Jan 28, 2018 at 11:50

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: So you'd want to find the variance of $\hat{a}+\epsilon$ which is estimated to be
nlm["CovarianceMatrix"][[1, 1]] + nlm["EstimatedVariance"]