# How to calculate margin of error in linear regression?

As a basic question in Regression Analysis, I wanted to ask how can I calculate Margin of Error when I fit a straight line to a set of data. Assume that I have variation of parameter $A$ as a function of parameter $B$, then $A=mB \pm e$, where $m$ is the tangent of the fitted straight line and $e$ is what i'm looking for!

Part of my data:

{{39.7678, 2320.3}, {30.8438, 1614.21}, {125.846, 3078.81}, {55.2345, \ 1947.98}, {22.0671, 972.995}, {30.1827, 701.99}, {29.5734, 837.784}, \ {24.6913, 1134.23}, {27.2493, 918.887}, {62.7684, 4535.07}, {101.449, \ 5499.83}, {125.248, 6513.04}, {187.409, 6257.72}, {174.138, 5243.63}, \ {120.747, 3768.02}, {84.178, 3453.12}, {60.2404, 3075.15}, {63.8622, \ 3517.73}, {101.9, 7240.11}, {90.6265, 5706.74}, {100.897, 7353.84}, \ {159.316, 9867.36}, {109.798, 11471.2}, {104.311, 6924.54}, {82.7057, \ 6339.06}, {140.205, 6555.52}, {173.469, 8644.27}, {138.432, 9655.86}, \ {95.2955, 5643.33}, {64.563, 3848.77}, {50.7936, 4733.24}, {34.776, 0. - 2707.89 I}, {25.3775, 6158.}}

-

You are asking to estimate the typical size of the error. This is usually done by estimating the variance of the random variable $e$ and taking its square root. It is often called the "root mean square error".

These data look like they are formatted for Mathematica. A Mathematica 8 solution is

{a, b} = {{39.7678, 2320.3},...<31 values omitted>...{25.3775, 6158.}} // Transpose;
model = LinearModelFit[{{b} // Transpose, a}];
model["EstimatedVariance"] // Sqrt


This assumes the data are given in the sequence (a,b) and not (b,a)!

It is important to note that this model omits the usual constant term. If you intended that to be there, include it by inserting a column of 1's in the "design matrix":

m = {ConstantArray[1, Length[a]], b} // Transpose;
model = LinearModelFit[{m, a}];
model["EstimatedVariance"] // Sqrt


In either case the output is a single number: your "margin of error."

-
 I found a link in MathWorld, about margin of error: mathworld.wolfram.com/MarginofError.html – user3341 Feb 23 '11 at 0:22 @K-1 Good find. Note the "additional qualification" that is needed: there is no consensus on what a "margin of error" really is. Notice, too, that the Wikipedia article offers three definitions: en.wikipedia.org/wiki/Margin_of_error#Definition . Regardless, in simple situations (like your least squares fit) the "margin of error" ultimately is computed from an estimate of the sampling variance, as shown in my sample code. Its square root, the sampling standard deviation, is often multiplied by some constant. – whuber♦ Feb 23 '11 at 0:25 1- By the way, it looks that the calculated error margin has no relation with "Confidence Level". For example if you add the option ConfidenceInterval->0.9 to the LinearModelFit there would be no change in the root mean square error. – user3341 Feb 23 '11 at 1:43 2- Instead of using LinearModelFit to calculate the constant in the fit, I suggest to use NonlinearModelFit. It enables you to define the pattern of fit. – user3341 Feb 23 '11 at 1:49 @K-1 (1) It all depends on what you are looking for confidence in. This is yet another reason why "margin of error" is so vague. You are better off computing an actual confidence interval (for a parameter or a prediction) and calling it that. (2) NonlinearModelFit is not what you specified in your question: you explicitly wrote down a linear model. – whuber♦ Feb 23 '11 at 14:48

It seems that you want the residuals of linear regression without the intercept term where dependent variable is $A$ and the independent variable is $B$. This can be done with various statistical packages. Here is the implementation in R.

aa <-"(39.7678, 2320.3}, {30.8438, 1614.21}, {125.846, 3078.81}, {55.2345,  1947.98}, {22.0671, 972.995}, {30.1827, 701.99}, {29.5734, 837.784},  {24.6913, 1134.23}, {27.2493, 918.887}, {62.7684, 4535.07}, {101.449,  5499.83}, {125.248, 6513.04}, {187.409, 6257.72}, {174.138, 5243.63},  {120.747, 3768.02}, {84.178, 3453.12}, {60.2404, 3075.15}, {63.8622,  3517.73}, {101.9, 7240.11}, {90.6265, 5706.74}, {100.897, 7353.84},  {159.316, 9867.36}, {109.798, 11471.2}, {104.311, 6924.54}, {82.7057,  6339.06}, {140.205, 6555.52}, {173.469, 8644.27}, {138.432, 9655.86},  {95.2955, 5643.33}, {64.563, 3848.77}, {50.7936, 4733.24}, {34.776, 2707.89 }, {25.3775, 6158)"

aa <- gsub(" *","",aa)

aa <- gsub("[}],[{]",");(",aa)

df<-t(sapply(strsplit(aa,";")[[1]],function(l)eval(parse(text=paste("c",l,sep="")))))

rownames(df) <-NULL

colnames(df)<-c("A","B")


All the code above is to read in your data into R. It would helped a lot if you simply provided the link to txt or csv file. Note that I fixed some errors in the last two sets of data.

Here are the residuals:

 residuals(lsfit(y=df[,"A"],x=df[,"B"],intercept=FALSE))

[1]   1.4348447   4.1759382  74.9819367  23.0525277   5.9925519
[6]  18.5853413  15.7326309   5.9530386  12.0686534 -12.1540823
[11]  10.5880176  17.6480900  84.0271524  87.5096181  58.4967124
[16]  27.1300775   9.4368024   5.7468757 -17.7116076  -3.6527783
[21] -20.5935056  -3.6995884 -79.7141307 -10.0871737 -22.0199405
[26]  31.9032911  30.6597014 -21.0894627   2.0637986   0.9786668
[31] -27.4027873  -9.9602093 -76.3569045


Note that these residuals have non-zero mean.

res <- residuals(lsfit(y=df[,"A"],x=df[,"B"],intercept=FALSE))
mean(res)

[1] 6.779518


So clearly you need to include intercept, since otherwise you cannot write $\pm e$. As @whuber pointed out, you probably need root mean square error. For the model with the intercept you can calculate is as follows:

> res <- residuals(lsfit(y=df[,"A"],x=df[,"B"],intercept=TRUE))
> res
[1] -15.716076 -15.933997  61.009708   4.341324 -16.804527  -5.347439
[7]  -7.631077 -16.168352 -10.955176 -20.023554   6.761571  18.067708
[13]  83.376799  82.609513  47.412760  14.726472  -4.550764  -6.385968
[19] -14.245054  -6.612126 -16.650343  10.776996 -58.516325  -7.943081
[25] -22.329422  32.500930  40.010676  -7.499212  -1.161282  -9.766886
[31] -34.441787 -25.486854 -77.425155

> sqrt(sum(res^2)/(length(res)-2))
[1] 34.70359
>


Note that you can use functions lm and summary to get all the proper statistics of linear regression:

> summary(lm(A~B,data=data.frame(df)))

Call:
lm(formula = A ~ B, data = data.frame(df))

Residuals:
Min      1Q  Median      3Q     Max
-77.425 -16.168  -7.499  10.777  83.377

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 26.874614  12.027882   2.234   0.0328 *
B            0.012330   0.002169   5.685 3.02e-06 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 34.7 on 31 degrees of freedom
Multiple R-squared: 0.5104, Adjusted R-squared: 0.4946
F-statistic: 32.32 on 1 and 31 DF,  p-value: 3.018e-06

-