I'm trying to calculate the uncertainty of the coefficients in a linear fit, but I'm not sure how to do it. on the one hand, I can extract the error with matlab statistical tools but on the other hand I need to take into account the uncertainty of my measurement (Let's say I measure pressure in a chamber and I have some degree of uncertainty that based on the measurement tool), and I dont know how to put it in my statistical calculations in matlab.

the code that I'm using now to evaluate the uncertainties is: [z,s]=polyfit(t,lnP,1); ste =sqrt(diag(inv(s.R)*inv(s.R')).*s.normr.^2./s.df);

Thank you for your help!

  • $\begingroup$ Can't you use the confidence interval for this purpose? $\endgroup$ – Kyle Kanos Jan 11 '17 at 15:32
  • $\begingroup$ Would Cross Validated be a better home for this question? $\endgroup$ – Qmechanic Jan 11 '17 at 15:47
  • 2
    $\begingroup$ I'm voting to close this question as off-topic because it is about programming/coding. $\endgroup$ – sammy gerbil Jan 12 '17 at 2:46
  • $\begingroup$ I think this question could be on topic if it is made to be more about the uncertainty in linear fits due to measurement errors (which I still think the CI can be used here) than 'how to do it in Matlab.' But at the same time, I think Qmechanics might be right and that this question might be better suited for Cross Validated. $\endgroup$ – Kyle Kanos Jan 12 '17 at 11:07
  • $\begingroup$ Hi! Yes, this is what I mean - Kyle understood my question when he commented : "How are uncertainties in the of measurement accounted for in linear fits" . I'm measuring the pressure inside a chamber that connects to a rotary pump. from the linear fit of ln(P) I can calculate the effective pump speed. My problem is how to evaluate the uncertainty of this result. On the one hand I can calculate it from the linear fit (Like Kevin estimated in his comment) $\endgroup$ – Noam Chai Jan 12 '17 at 11:26

let's assume you have N data points of the form $(x_i,y_i\pm s_i)$, and you want some function of x to approximate y. I'll assume a more general form for the fitting function, but it will be applicable (and much simpler) in the case of a linear fit. $$f(x)=\sum_{i=1}^n c_i x^{e_i}$$ For a linear fit, the exponents ${e_i}$ will just be zero and one, and you'll need to optimize the coefficients ${c_i}$.

I assume you're doing a least-squares fit, so the quantity you'll want to minimize is $\chi^2$: $$\chi^2=\sum_{i=1}^N\left(\frac{y_i-f(x_i)}{s_i}\right)^2$$

To minimize this with respect to $c_i$, you'll need to find a stationary point such that $\frac{\partial\chi^2}{\partial c_i}=0$ for all $c_i$. This can be cast as a system of linear equations: $$\frac{\partial\chi^2}{\partial c_k}=2\sum_{i=1}^N\left(\frac{y_i-f(x_i)}{s_i}\right)\left(-\frac{x_i^{e_k}}{s_i}\right)=0$$ $$\sum_{i=1}^N\left(\frac{y_i-\sum_{j=1}^n c_j x_i^{e_j}}{s_i}\right)\left(-\frac{x_i^{e_k}}{s_i}\right)=0$$ $$\sum_{i=1}^N\left(\frac{x_i^{e_k}\sum_{j=1}^n c_j x_i^{e_j}}{s_i^2}\right)=\sum_{i=1}^N\left(\frac{y_i x_i^{e_k}}{s_i^2}\right)$$ $$\sum_{i=1}^N\sum_{j=1}^n\left(\frac{x_i^{e_k+e_j}}{s_i^2}\right)c_j=\sum_{i=1}^N\left(\frac{y_i x_i^{e_k}}{s_i^2}\right)$$ $$\sum_{j=1}^n A_{kj}c_j=b_k$$ $$c_k=\sum_{j=1}^n A_{kj}^{-1}b_j$$ Where $$A_{kj}=\sum_{i=1}^N\left(\frac{x_i^{e_k+e_j}}{s_i^2}\right)$$ $$b_k=\sum_{i=1}^N\left(\frac{y_i x_i^{e_k}}{s_i^2}\right)$$

The variance $\sigma_i^2$ of each coefficient $c_i$ can be estimated by (see equation 9) \begin{align*} \sigma_i^2 &= \sum_{j=1}^N s_j^2 \left(\frac{\partial c_i}{\partial y_j}\right)^2\\ &= \sum_{j=1}^N s_j^2 \left(\sum_{k=1}^n\frac{\partial A_{ik}^{-1}b_k}{\partial y_j}\right)^2\\ &= \sum_{j=1}^N s_j^2 \left(\sum_{k=1}^n A_{ik}^{-1}\frac{\partial b_k}{\partial y_j}\right)^2\\ &= \sum_{j=1}^N s_j^2 \left(\sum_{k=1}^n A_{ik}^{-1}\frac{x_j^{e_k}}{s_j^2}\right)^2 \end{align*}

So as long as $A$ is invertible, you can use the procedure above to find uncertainties in your coefficients.

the $\sigma_i^2$ actually simplify a bit further:

\begin{align*} \sigma_i^2 &= \sum_{j=1}^N s_j^2 \left(\sum_{k=1}^n A_{ik}^{-1}\frac{x_j^{e_k}}{s_j^2}\right)^2\\ &=\sum_{j=1}^N s_j^2 \left(\sum_{k=1}^n A_{ik}^{-1}\frac{x_j^{e_k}}{s_j^2}\right)\left(\sum_{l=1}^n A_{il}^{-1}\frac{x_j^{e_l}}{s_j^2}\right)\\ &=\sum_{k=1}^n\sum_{l=1}^n A_{ik}^{-1} A_{il}^{-1}\sum_{j=1}^N\frac{x_j^{e_l+e_k}}{s_j^2}\\ &=\sum_{k=1}^n\sum_{l=1}^n A_{ik}^{-1} A_{il}^{-1}A_{lk}\\ &=\sum_{k=1}^n A_{ik}^{-1} \delta_{ik}\\ \sigma_i^2&=A_{ii}^{-1} \end{align*} So your function of x is now $$f(x)=\sum_{i=1}^n (c_i \pm \sigma_i)x^{e_i}$$ The uncertainty in the evaluation of the function at some point $x_j$ can be estimated similarly to how the uncertainty in coefficients was calculated.

\begin{align*} \sigma_{f(x_j)}^2 &=\sum_{i=1}^n\left(\frac{\partial f(x_j)}{\partial c_i}\right)^2 \sigma_i^2\\ &=\sum_{i=1}^n\left(x_j^{e_i}\right)^2 \sigma_i^2 \end{align*}

Given some point $x_i$, the value given by the fit will be $f(x_i)$, and the uncertainty (one std. deviation) will be $\sigma_{f(x_i)}$ (as defined above)

  • $\begingroup$ This seems to answer the question "How does one do a linear fit?" which isn't what OP asked. They asked "How are uncertainties in the of measurement accounted for in linear fits?" $\endgroup$ – Kyle Kanos Jan 12 '17 at 11:01
  • $\begingroup$ @Kyle Kanos He asked how to calculate the uncertainties in the coefficients of a linear fit, which I explained how to do. I included the steps for the rest of the fitting process because I thought the expressions for the uncertainties would make more sense in that context. $\endgroup$ – kevin Jan 12 '17 at 13:50
  • $\begingroup$ That is what the first sentence says, but OP has clarified that my interpretation of the question is what he was aiming for: he has error in data-taking and wants that accounted for in his statistical analysis. $\endgroup$ – Kyle Kanos Jan 12 '17 at 13:56
  • $\begingroup$ I don't have enough reputation to comment on the original question, but I don't see what the OP is asking for. Given uncertainties in the $y$ values (error in data-taking), I show how to determine uncertainties in the coefficients and subsequently how to determine the uncertainty in the value of the linear fit at any point (error in model). The $\chi^2$ value can also be used in some kind of chi-squared test (error in analysis?). $\endgroup$ – kevin Jan 12 '17 at 22:13

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