0
$\begingroup$

I have a following data which follows a linear model with intercept zero, i.e. $$y = mx $$

----------------------
| x-value | y-values |
|--------------------|
|  1      | 0.9±0.2  |
|  1.5    | 1.2±0.2  |
|  2      | 1.5±0.2  |
|  2.5    | 2.0±0.2  |
|  3      | 2.6±0.2  |
----------------------

How do I find the maximum and minimum value for $m$.

$\endgroup$
7
$\begingroup$

Looking at help page for LinearModelFit,

Using VarianceEstimatorFunction->(1&) and Weights->{1/Δy_1^2,1/Δy_2^2,…}, Δy_i is treated as the known uncertainty of measurement y_i and parameter standard errors are effectively computed only from the weights.

data = {{1, .9}, {1.5, 1.2}, {2., 1.5}, {2.5, 2.0}, {3., 2.6}};
errors = {.2, .2, .2, .2, .2};
lmod = LinearModelFit[data, x, x, VarianceEstimatorFunction -> (1 &), 
  Weights -> 1/errors^2]

enter image description here

You can then get the errors for the slope and intercept via

lmod["ParameterTable"]

enter image description here

So the range of slope values is

.84 + .126491 {1, -1}
(* {0.966491, 0.713509} *)

You can look at the range of results using an ErrorListPlot

Needs["ErrorBarPlots`"]
Show[
 ErrorListPlot[Transpose[{data[[All, 1]], data[[All, 2]], errors}], 
  PlotStyle -> Red],
 Plot[{-.04 + .84 x,
   -.04 + .268328 + (.84 - .126491) x,
   -.04 - 0.268328 + (.84 + .126491) x},
  {x, 0, 5}]
 ]

enter image description here

$\endgroup$
  • $\begingroup$ thank you for answer. Could it be possible to add error values for $x$ values too? $\endgroup$ – Santosh Linkha Apr 14 '16 at 10:44
  • $\begingroup$ I don't think you can do that out of the box, but this post suggests taking the root sum square of the two errors (x and y) and use that as an effective y error. This post looks more rigorous, but also more complicated. $\endgroup$ – Jason B. Apr 14 '16 at 10:51
  • $\begingroup$ all right thanks. I was looking if there were easy method. $\endgroup$ – Santosh Linkha Apr 14 '16 at 10:55
3
$\begingroup$

Warning: The result at the end is very anticlimactic.

While you’ve stated that your model has an intercept of zero the model that @JasonB is assuming with the particular use of LinearModelFit is the following:

$$y_i=a+b x_i+ϵ_i$$

with $ϵ_i\sim N(0.2^2)$. Besides having an intercept (potentially not equal to zero), this model assumes that there is only measurement error (without any lack-of-fit error) and that the standard deviation of the measurement error is known (and therefore doesn’t have to be estimated).

If the plus-or-minus numbers that you give for the y-values are really known measurement standard deviations, possibly a more appropriate model is the following:

$$y_i=b x_i+ϵ_i+γ_i$$

where $ϵ_i\sim N(0.2^2)$ are the measurement errors and $γ_i\sim N(0,σ^2)$ are the lack-of-fit errors (i.e., deviations from the line when there is no measurement error) with the errors being independent of each other. (And I acknowledge that performing this more complicated model on only 5 data points is not really recommended. I offer this only as an example that would be used with many more data points.)

While in this case the estimate of the slope won’t change much at all if one fits

$$y_i=b x_i+ϵ_i+γ_i$$ or one fits $$y_i=b x_i+ϵ_i$$ It is the estimate of the standard error of the slope that has the potential to change. One is trying to characterize and use all the information that is available to obtain appropriate estimates of both the parameters and an associated measure of precision.

LinearModelFit and NonlinearModelFit only allow a single error term with a specific functional form so one way to fit this model with the structure of the model and the data available is to find the maximum likelihood estimates using the LogLikelihood function. The variance of $y_i |x_i$ ($y_i$ given $x_i$) is $0.2^2+σ^2$ so the log of the likelihood function is

logL = Simplify[LogLikelihood[NormalDistribution[0, ((1/5)^2 + σ^2)^(1/2)],
   data[[All, 2]] - b data[[All, 1]]]]
FindMaximum[logL && σ >= 0 , {{b, 0.84}, {σ, 0}}]
(* {2.7525, {b -> 0.84, σ -> 0.}} *)

So our estimate of $\sigma$ is zero. To look into this further suppose we ignored the “known” measurement standard deviation and use the standard approach for fitting a line through the origin:

lm = LinearModelFit[data, x, x, IncludeConstantBasis -> False];
lm["EstimatedVariance"]^0.5
(* 0.110554 *)
lm["ParameterTableEntries"]
(* {{0.822222,0.0233069,35.27811,3.853070 * 10^-6}} *)
lm["ParameterConfidenceIntervals"]
(* {{0.757512,0.886932}} *)

We get an estimate for the the slope of 0.822222 and an estimated standard deviation about the line of 0.110554 which is about half of the “known” standard deviation. So you might not know what you think you know. (If what you presented in your original question was +/- 2 standard deviations, then you're right on the money.)

The estimated standard error for the slope is 0.0233069. A 95% confidence interval for the slope is {0.757512, 0.886932}. If want just to use the estimate plus-or-minus one standard error, then that should be labeled with very similar words. (There is no such thing as a "min-max for slope".)

$\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.