|
From: <pl...@pi...> - 2011-05-21 08:15:37
|
Hi, 'help fit' reports that the fit command uses Levenberg–Marquardt algo to do the fit. I think this raises an important question that very often over-looked by many users of least-squares techniques even at maths PhD level. Such techniques often only optimised the y error rather than the perpendicular error from the line. This is implicitly assuming y uncertainty >> x uncertainty. While this condition is often satisfied in a controlled experiment there are many situations where this is not applicable and gets totally overlooked. A common case is scatter plots which are frequently used to seek a relations between two quantities , each with significant errors / uncertainties. In this situation the fitted line is "wrong". In fact it's the application that is wrong , hence the wrong result. This may or may not be apparent to the eye. I have seen this happen so many times (including once in a PhD thesis report!) that I think it needs a serious health warning in the doc. "Warning: using least-squares inappropriately can seriously damage your reputation". ;) Firstly , could you confirm the basis on which this algo is applied in gnuplot? Does it only optimise vertical y residuals? Thanks, Peter. |
|
From: Daniel J S. <dan...@ie...> - 2011-05-21 21:04:51
|
On 05/21/2011 03:15 AM, pl...@pi... wrote: > Hi, > > 'help fit' reports that the fit command uses Levenberg–Marquardt algo to > do the fit. > > I think this raises an important question that very often over-looked by > many users of least-squares techniques even at maths PhD level. > > Such techniques often only optimised the y error rather than the > perpendicular error from the line. This is implicitly assuming y > uncertainty>> x uncertainty. While this condition is often satisfied > in a controlled experiment there are many situations where this is not > applicable and gets totally overlooked. > > A common case is scatter plots which are frequently used to seek a > relations between two quantities , each with significant errors / > uncertainties. > > In this situation the fitted line is "wrong". In fact it's the > application that is wrong , hence the wrong result. This may or may not > be apparent to the eye. > > I have seen this happen so many times (including once in a PhD thesis > report!) that I think it needs a serious health warning in the doc. > > "Warning: using least-squares inappropriately can seriously damage your > reputation". ;) > > Firstly , could you confirm the basis on which this algo is applied in > gnuplot? Does it only optimise vertical y residuals? Often it is an assumption that the independent variables are exact measurements. Not true, typically, but if the variance is small and homoscedastic, the two can probably be lumped together. I.e., we are searching for a relationship: Y = f(X + eps1) + eps2 ~= f(X) + C eps1 + residual + eps2 ~= f(X) + (C eps1 + eps2) where hopefully the residual due to nonlinearity of the relationship is small compared to other randomness. It's up to the user's judgment and knowledge of the application to determine that. Anyway, your point is true of most software packages: details are so often lacking. That's why it would be nice to have a set of white papers to go along with the software so that people know exactly what the algorithm is, both for the benefit of the user and other developers. Most of the time it is "here's a hunk of code, use it at your own risk". Dan |
|
From: <pl...@pi...> - 2011-05-22 18:14:26
|
On 05/21/11 23:04, Daniel J Sebald wrote: > On 05/21/2011 03:15 AM, pl...@pi... wrote: >> Hi, >> >> 'help fit' reports that the fit command uses Levenberg–Marquardt algo to >> do the fit. >> >> I think this raises an important question that very often over-looked by >> many users of least-squares techniques even at maths PhD level. >> >> Such techniques often only optimised the y error rather than the >> perpendicular error from the line. This is implicitly assuming y >> uncertainty>> x uncertainty. While this condition is often satisfied >> in a controlled experiment there are many situations where this is not >> applicable and gets totally overlooked. >> >> A common case is scatter plots which are frequently used to seek a >> relations between two quantities , each with significant errors / >> uncertainties. >> >> In this situation the fitted line is "wrong". In fact it's the >> application that is wrong , hence the wrong result. This may or may not >> be apparent to the eye. >> >> I have seen this happen so many times (including once in a PhD thesis >> report!) that I think it needs a serious health warning in the doc. >> >> "Warning: using least-squares inappropriately can seriously damage your >> reputation". ;) >> >> Firstly , could you confirm the basis on which this algo is applied in >> gnuplot? Does it only optimise vertical y residuals? > > Often it is an assumption that the independent variables are exact > measurements. Not true, typically, but if the variance is small and > homoscedastic, the two can probably be lumped together. I.e., we are > searching for a relationship: > > Y = f(X + eps1) + eps2 > ~= f(X) + C eps1 + residual + eps2 > ~= f(X) + (C eps1 + eps2) > > where hopefully the residual due to nonlinearity of the relationship is > small compared to other randomness. It's up to the user's judgment and > knowledge of the application to determine that. > > Anyway, your point is true of most software packages: details are so > often lacking. That's why it would be nice to have a set of white > papers to go along with the software so that people know exactly what > the algorithm is, both for the benefit of the user and other developers. > Most of the time it is "here's a hunk of code, use it at your own risk". > > Dan > Thanks, I've read up on that algo and clearly this is only doing NLLS on y residuals. So my suggestion is that this is made abundantly clear in the help text. I'm not suggesting your "while paper" but just some comment to the effect that 'fit' will not give correct results if there are non negligible errors in x values. It absolutely amazes me how few people realise this , even highly qualified ones, so this is not some pedantic nicety. Most people seem to think once they've heard of doing a least squares fit that's all there is to it and it's some magical formula that works for all cases. I suggest modifying the first paragraph of help fit with something like the following: >> The `fit` command can fit a user-supplied expression to a set of data points (x,z) or (x,y,z), using an implementation of the nonlinear least-squares (NLLS) Marquardt-Levenberg algorithm. Any user-defined variable occurring in the expression may serve as a fit parameter, but the return type of the expression must be real. >> new >> The `fit` command can fit a user-supplied expression to a set of data points (x,z) or (x,y,z), using an implementation of the non-linear least-squares (NLLS) Marquardt-Levenberg algorithm. This algorithm optimises y residuals only and carries the implicit assumption that error/uncertainties in x are negligible. If that is not the case the fit may succeed but will give wrong results. Any user-defined variable occurring in the expression may serve as a fit parameter, but the return type of the expression must be real. >> regards, Peter. |
|
From: Daniel J S. <dan...@ie...> - 2011-05-23 00:55:19
|
On 05/22/2011 01:14 PM, pl...@pi... wrote: > On 05/21/11 23:04, Daniel J Sebald wrote: >> On 05/21/2011 03:15 AM, pl...@pi... wrote: >>> Hi, >>> >>> 'help fit' reports that the fit command uses Levenberg–Marquardt algo to >>> do the fit. >>> >>> I think this raises an important question that very often over-looked by >>> many users of least-squares techniques even at maths PhD level. >>> >>> Such techniques often only optimised the y error rather than the >>> perpendicular error from the line. This is implicitly assuming y >>> uncertainty>> x uncertainty. While this condition is often satisfied >>> in a controlled experiment there are many situations where this is not >>> applicable and gets totally overlooked. >>> >>> A common case is scatter plots which are frequently used to seek a >>> relations between two quantities , each with significant errors / >>> uncertainties. >>> >>> In this situation the fitted line is "wrong". In fact it's the >>> application that is wrong , hence the wrong result. This may or may not >>> be apparent to the eye. >>> >>> I have seen this happen so many times (including once in a PhD thesis >>> report!) that I think it needs a serious health warning in the doc. >>> >>> "Warning: using least-squares inappropriately can seriously damage your >>> reputation". ;) >>> >>> Firstly , could you confirm the basis on which this algo is applied in >>> gnuplot? Does it only optimise vertical y residuals? >> >> Often it is an assumption that the independent variables are exact >> measurements. Not true, typically, but if the variance is small and >> homoscedastic, the two can probably be lumped together. I.e., we are >> searching for a relationship: >> >> Y = f(X + eps1) + eps2 >> ~= f(X) + C eps1 + residual + eps2 >> ~= f(X) + (C eps1 + eps2) >> >> where hopefully the residual due to nonlinearity of the relationship is >> small compared to other randomness. It's up to the user's judgment and >> knowledge of the application to determine that. >> >> Anyway, your point is true of most software packages: details are so >> often lacking. That's why it would be nice to have a set of white >> papers to go along with the software so that people know exactly what >> the algorithm is, both for the benefit of the user and other developers. >> Most of the time it is "here's a hunk of code, use it at your own risk". >> >> Dan >> > > Thanks, I've read up on that algo and clearly this is only doing NLLS on > y residuals. > > So my suggestion is that this is made abundantly clear in the help text. > I'm not suggesting your "while paper" but just some comment to the > effect that 'fit' will not give correct results if there are non > negligible errors in x values. > > It absolutely amazes me how few people realise this , even highly > qualified ones, so this is not some pedantic nicety. > > Most people seem to think once they've heard of doing a least squares > fit that's all there is to it and it's some magical formula that works > for all cases. > > I suggest modifying the first paragraph of help fit with something like > the following: > > >> > The `fit` command can fit a user-supplied expression to a set of data > points > (x,z) or (x,y,z), using an implementation of the nonlinear least-squares > (NLLS) Marquardt-Levenberg algorithm. Any user-defined variable > occurring in > the expression may serve as a fit parameter, but the return type of the > expression must be real. > >> > > new > >> > The `fit` command can fit a user-supplied expression to a set of data > points > (x,z) or (x,y,z), using an implementation of the non-linear least-squares > (NLLS) Marquardt-Levenberg algorithm. This algorithm optimises y > residuals only and > carries the implicit assumption that error/uncertainties in x are > negligible. If that is not the case the fit may succeed but will give > wrong results. I'm OK with the change, but for a few things. First, you stated 'optimizes y' and 'uncertainties in x', but please check that this precisely describes the algorithm. The beginning of the paragraph lists (x, z) or (x, y, z). The first expression has no 'y', so what is optimized in that case? The second expression has (x, y, z), so is it only the 'x' in that case that is assumed exact? Or is it both 'x' and 'y'? Perhaps that sentence should be, "This algorithm optimises z residuals only and carries the implicit assumption that error/uncertainties in x (and y) are negligible." Second, I would hold off using the statement "wrong results" and maybe use "a poor fit". Saying the result is wrong means the algorithm is broken, but it just provides numbers. The user is using the wrong tool. Also, "wrong" means a bad fit in this context and one can get a bad fit and misinterpret results even if uncertainty in x is small. Third, if there is an alternative approach for the case when x uncertainty can't be ignored, keep the discussion going and maybe we can add such a feature to the list of items we'd like to add in the future. Dan |
|
From: <pl...@pi...> - 2011-05-23 08:39:02
|
On 05/23/11 02:55, Daniel J Sebald wrote: > On 05/22/2011 01:14 PM, pl...@pi... wrote: >> On 05/21/11 23:04, Daniel J Sebald wrote: >>> On 05/21/2011 03:15 AM, pl...@pi... wrote: >>>> Hi, >>>> >>>> 'help fit' reports that the fit command uses Levenberg–Marquardt algo to >>>> do the fit. >>>> >>>> I think this raises an important question that very often over-looked by >>>> many users of least-squares techniques even at maths PhD level. >>>> >>>> Such techniques often only optimised the y error rather than the >>>> perpendicular error from the line. This is implicitly assuming y >>>> uncertainty>> x uncertainty. While this condition is often satisfied >>>> in a controlled experiment there are many situations where this is not >>>> applicable and gets totally overlooked. >>>> >>>> A common case is scatter plots which are frequently used to seek a >>>> relations between two quantities , each with significant errors / >>>> uncertainties. >>>> >>>> In this situation the fitted line is "wrong". In fact it's the >>>> application that is wrong , hence the wrong result. This may or may not >>>> be apparent to the eye. >>>> >>>> I have seen this happen so many times (including once in a PhD thesis >>>> report!) that I think it needs a serious health warning in the doc. >>>> >>>> "Warning: using least-squares inappropriately can seriously damage your >>>> reputation". ;) >>>> >>>> Firstly , could you confirm the basis on which this algo is applied in >>>> gnuplot? Does it only optimise vertical y residuals? >>> >>> Often it is an assumption that the independent variables are exact >>> measurements. Not true, typically, but if the variance is small and >>> homoscedastic, the two can probably be lumped together. I.e., we are >>> searching for a relationship: >>> >>> Y = f(X + eps1) + eps2 >>> ~= f(X) + C eps1 + residual + eps2 >>> ~= f(X) + (C eps1 + eps2) >>> >>> where hopefully the residual due to nonlinearity of the relationship is >>> small compared to other randomness. It's up to the user's judgment and >>> knowledge of the application to determine that. >>> >>> Anyway, your point is true of most software packages: details are so >>> often lacking. That's why it would be nice to have a set of white >>> papers to go along with the software so that people know exactly what >>> the algorithm is, both for the benefit of the user and other developers. >>> Most of the time it is "here's a hunk of code, use it at your own risk". >>> >>> Dan >>> >> >> Thanks, I've read up on that algo and clearly this is only doing NLLS on >> y residuals. >> >> So my suggestion is that this is made abundantly clear in the help text. >> I'm not suggesting your "while paper" but just some comment to the >> effect that 'fit' will not give correct results if there are non >> negligible errors in x values. >> >> It absolutely amazes me how few people realise this , even highly >> qualified ones, so this is not some pedantic nicety. >> >> Most people seem to think once they've heard of doing a least squares >> fit that's all there is to it and it's some magical formula that works >> for all cases. >> >> I suggest modifying the first paragraph of help fit with something like >> the following: >> >> >> >> The `fit` command can fit a user-supplied expression to a set of data >> points >> (x,z) or (x,y,z), using an implementation of the nonlinear least-squares >> (NLLS) Marquardt-Levenberg algorithm. Any user-defined variable >> occurring in >> the expression may serve as a fit parameter, but the return type of the >> expression must be real. >> >> >> >> new >> >> >> The `fit` command can fit a user-supplied expression to a set of data >> points >> (x,z) or (x,y,z), using an implementation of the non-linear least-squares >> (NLLS) Marquardt-Levenberg algorithm. This algorithm optimises y >> residuals only and >> carries the implicit assumption that error/uncertainties in x are >> negligible. If that is not the case the fit may succeed but will give >> wrong results. > > I'm OK with the change, but for a few things. > > First, you stated 'optimizes y' and 'uncertainties in x', but please > check that this precisely describes the algorithm. The beginning of the > paragraph lists (x, z) or (x, y, z). The first expression has no 'y', > so what is optimized in that case? The second expression has (x, y, z), > so is it only the 'x' in that case that is assumed exact? Or is it both > 'x' and 'y'? Perhaps that sentence should be, "This algorithm optimises > z residuals only and carries the implicit assumption that > error/uncertainties in x (and y) are negligible." very good point , I was trying to keep it brief (and perhaps dumb it down) but it should correctly refer to dependent and independent variables. I was concerned that may mean the warning was lost on many users. Perhaps "independant variable (eg. x axis)" or similar would be better. > > Second, I would hold off using the statement "wrong results" and maybe > use "a poor fit". Saying the result is wrong means the algorithm is > broken, but it just provides numbers. The user is using the wrong tool. > Also, "wrong" means a bad fit in this context and one can get a bad > fit and misinterpret results even if uncertainty in x is small. No, I think wrong result is exactly the case. Wrong result does not mean the algo is "broken" it means the wrong tool was used and the result obtained was wrong. The result is not "poor" it is wrong. Watering down the language does not shift the blame. Maybe some wording like "wrong technique" could help underline the cause of the problem but I think "wrong" definitely needs to be in there. There's no hedging around the fact, it can be considerably wrong. > > Third, if there is an alternative approach for the case when x > uncertainty can't be ignored, keep the discussion going and maybe we can > add such a feature to the list of items we'd like to add in the future. > > Dan > Yes , I would love to propose that. Some kind of total least squares may be an option I don't know enough about that to make a concrete proposal. I suspect it may be too variable to include as a turn key option like fit but I would really like that option if it is possible. thanks for your thoughtful comments. |
|
From: Ethan M. <merritt@u.washington.edu> - 2011-05-24 00:08:17
|
On Monday, May 23, 2011 01:12:23 am pl...@pi... wrote: > > if there is an alternative approach for the case when x > > uncertainty can't be ignored, keep the discussion going and maybe we can > > add such a feature to the list of items we'd like to add in the future. > > > > Dan > > Yes , I would love to propose that. Some kind of total least squares may > be an option I don't know enough about that to make a concrete proposal. It is standard in my field to use instead a maximum likelihood residual that allows for separate error distributions on x and y. The maximum likelihood treatment reduces to a weighted least squares treatment if and only if the distribution of errors on both x and y are (1) Gaussian and (2) of equal magnitude. Gnuplot allows for input of precalculated non-uniform weights on y, which partially addresses (2). But there is no provision for non-Gaussian errors on y, and no provision for any error model at all on x. If the errors in your data do not follow the simple Gaussian model, then almost certainly you would be better off using maximum likelihood rather than least-squares. But in order to do so you need first a model for the error distributions. That may be obvious for any particular experiment or source of data, but it's difficult to impossible for a general-purpose program to figure it out for you. You're the one who knows the source of the data, so it's up to you to provide an appropriate explicit error model along with the data. The point is that switching to a better minimization residual is more than just a matter of changing the internal code. It would require a more complex description of your data that includes error models for both the independent and dependent variables. NB: When I say "x", I really mean the full set of independent variables [x1,x2,x3,...]. "y" is the single dependent variable estimated by f(x1,x2,x3,...) Ethan -- Ethan A Merritt Biomolecular Structure Center, K-428 Health Sciences Bldg University of Washington, Seattle 98195-7742 |
|
From: <pl...@pi...> - 2011-05-24 06:53:28
|
On 05/24/11 02:08, Ethan Merritt wrote: > On Monday, May 23, 2011 01:12:23 am pl...@pi... wrote: >>> if there is an alternative approach for the case when x >>> uncertainty can't be ignored, keep the discussion going and maybe we can >>> add such a feature to the list of items we'd like to add in the future. >>> >>> Dan >> >> Yes , I would love to propose that. Some kind of total least squares may >> be an option I don't know enough about that to make a concrete proposal. > > It is standard in my field to use instead a maximum likelihood residual > that allows for separate error distributions on x and y. The maximum > likelihood treatment reduces to a weighted least squares treatment > if and only if the distribution of errors on both x and y are > (1) Gaussian and (2) of equal magnitude. > > Gnuplot allows for input of precalculated non-uniform weights on y, > which partially addresses (2). But there is no provision for non-Gaussian > errors on y, and no provision for any error model at all on x. > > If the errors in your data do not follow the simple Gaussian model, > then almost certainly you would be better off using maximum likelihood > rather than least-squares. But in order to do so you need first a model > for the error distributions. That may be obvious for any particular > experiment or source of data, but it's difficult to impossible for a > general-purpose program to figure it out for you. You're the one who > knows the source of the data, so it's up to you to provide an appropriate > explicit error model along with the data. > > The point is that switching to a better minimization residual is more > than just a matter of changing the internal code. It would require a > more complex description of your data that includes error models for > both the independent and dependent variables. > > NB: When I say "x", I really mean the full set of independent variables > [x1,x2,x3,...]. "y" is the single dependent variable estimated by > f(x1,x2,x3,...) > > Ethan > Thanks for all that detail Ethan. It seems likely that this may be the job for external processing by some stats capable package rather than a job for a plotting tool. That brings us back to the original point of this thread : the idea of putting a specific warning about the limitations of least squares techniques used by gnuplot into the help text for "fit". I won't bore everyone with anecdotes but I never cease to be amazed by the widespread ignorance of this issue , even at PhD level. For the cost of a couple of lines it would be good if at least gnuplot user were made aware of it. best regards.Peter. |