Uncertainty analysis library?

newer
wrong backquote in haskell 2010...

Edward Amsden

20 Mar 2011 20 Mar '11

4:27 p.m.

Hi cafe, I'm looking for a library that provides an instance of Num, Fractional, Floating, etc, but carries uncertainty values through calculations. A scan of hackage didn't turn anything up. Does anyone know of a library like this? Thanks! -- Edward Amsden Student Computer Science Rochester Institute of Technology www.edwardamsden.com

Show replies by date

Henning Thielemann

20 Mar 20 Mar

4:34 p.m.

On Sun, 20 Mar 2011, Edward Amsden wrote:

...

Hi cafe,

I'm looking for a library that provides an instance of Num, Fractional, Floating, etc, but carries uncertainty values through calculations. A scan of hackage didn't turn anything up. Does anyone know of a library like this?

Do you mean exact bounds, that is, interval arithmetic? Or do you mean rough estimates of the uncertainty - this might be related to automatic differentiation.

Edward Kmett

4:38 p.m.

I have a package for interval arithmetic in hackage http://hackage.haskell.org/package/intervals-0.2.0 However it does not currently properly adjust the floating point rounding mode so containment isn't perfect. However, we are actively working on fixing up the Haskell MPFR bindings, which will let us reliably set rounding modes, cleaning up the interval arithmetic library to be just a little bit more pedantic. Due to the way GHC interacts with GMP this is a disturbingly difficult process. I have an unreleased library for working with Taylor models that builds on top of that and my automatic differentiation library, but without working MPFR bindings, it isn't sufficiently accurate for me to comfortably release. -Edward On Sun, Mar 20, 2011 at 4:27 PM, Edward Amsden wrote:

...

Hi cafe,

I'm looking for a library that provides an instance of Num, Fractional, Floating, etc, but carries uncertainty values through calculations. A scan of hackage didn't turn anything up. Does anyone know of a library like this?

Thanks!

-- Edward Amsden Student Computer Science Rochester Institute of Technology www.edwardamsden.com

_______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe

Tom Nielsen

5:46 p.m.

Interval arithmetic is of course not the same as uncertainty, although computer scientists like to pretend that is the case. (and uncertainty estimates do not have the be "rough".) In general the propagation of errors depends on whether the errors are independent or not. The rules are given in Taylor: An introduction to Error analysis (1997). Interval artihmetic corresponds to the worst case of non-independent and non-random errors. In the case of independent of random errors, you get: data Approximately a = a :+/-: a instance Num a => Num (Approximately a) where (m1 :+/-: err1) + (m2 :+/-: err2) = (m1+m2) :+/-: (sqrt(err1^2+err2^2) (m1 :+/-: err1) - (m2 :+/-: err2) = (m1-m2) :+/-: (sqrt(err1^2+err2^2) (m1 :+/-: err1) * (m2 :+/-: err2) = (m1*m2) :+/-: (sqrt((err1/m1)^2+(err2/m2)^2) the general rule is if y = f xs where xs :: [Approximately a], i.e f :: [Approximately a] -> Approximately a the error term= sqrt $ sum $ map (^2) $ map (\(ym :+/-: yerr) -> partial-derivative-of-yerr-with-respect-to-partial-ym * yerr) xs You can verify these things by running your calculation through soem sort of randomness monad (monte-carlo or random-fu packages) Anyways, I ended up not going down this route this because probabilistic data analysis gives you the correct error estimate without propagating error terms. Tom PS if you're a scientist and your accuracy estimate is on the same order as your rounding error, your are doing pretty well :-) At least in my field... On Sun, Mar 20, 2011 at 8:38 PM, Edward Kmett wrote:

...

I have a package for interval arithmetic in hackage http://hackage.haskell.org/package/intervals-0.2.0 However it does not currently properly adjust the floating point rounding mode so containment isn't perfect. However, we are actively working on fixing up the Haskell MPFR bindings, which will let us reliably set rounding modes, cleaning up the interval arithmetic library to be just a little bit more pedantic. Due to the way GHC interacts with GMP this is a disturbingly difficult process. I have an unreleased library for working with Taylor models that builds on top of that and my automatic differentiation library, but without working MPFR bindings, it isn't sufficiently accurate for me to comfortably release. -Edward

On Sun, Mar 20, 2011 at 4:27 PM, Edward Amsden wrote:

...
Hi cafe,

I'm looking for a library that provides an instance of Num, Fractional, Floating, etc, but carries uncertainty values through calculations. A scan of hackage didn't turn anything up. Does anyone know of a library like this?

Thanks!

-- Edward Amsden Student Computer Science Rochester Institute of Technology www.edwardamsden.com

_______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe

_______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe

Edward Amsden

7:59 p.m.

I'm actually a CS undergrad in a physics lab class. I have permission from my professor to use computer programs for analysis of lab data. I need to do calculations on data with uncertainty, but uncertainty analysis on many formulae in physics is rather tedious. I was hoping for something with instances for Num, Fractional, Floating, etc. that would allow me to combine two uncertain values and get a new value with uncertainty. I've been working on writing one myself and I don't find the concept hard, but it's a lot of effort that I don't want to duplicate if it's been done already. On Sun, Mar 20, 2011 at 5:46 PM, Tom Nielsen wrote:

...

Interval arithmetic is of course not the same as uncertainty, although computer scientists like to pretend that is the case. (and uncertainty estimates do not have the be "rough".)

In general the propagation of errors depends on whether the errors are independent or not. The rules are given in Taylor: An introduction to Error analysis (1997). Interval artihmetic corresponds to the worst case of non-independent and non-random errors. In the case of independent of random errors, you get:

data Approximately a = a :+/-: a

instance Num a => Num (Approximately a) where (m1 :+/-: err1) + (m2 :+/-: err2) = (m1+m2) :+/-: (sqrt(err1^2+err2^2) (m1 :+/-: err1) - (m2 :+/-: err2) = (m1-m2) :+/-: (sqrt(err1^2+err2^2) (m1 :+/-: err1) * (m2 :+/-: err2) = (m1*m2) :+/-: (sqrt((err1/m1)^2+(err2/m2)^2)

the general rule is

if y = f xs where xs :: [Approximately a], i.e f :: [Approximately a] -> Approximately a

the error term= sqrt $ sum $ map (^2) $ map (\(ym :+/-: yerr) -> partial-derivative-of-yerr-with-respect-to-partial-ym * yerr) xs

You can verify these things by running your calculation through soem sort of randomness monad (monte-carlo or random-fu packages) Anyways, I ended up not going down this route this because probabilistic data analysis gives you the correct error estimate without propagating error terms.

Tom

PS if you're a scientist and your accuracy estimate is on the same order as your rounding error, your are doing pretty well :-) At least in my field...

On Sun, Mar 20, 2011 at 8:38 PM, Edward Kmett wrote:

...
I have a package for interval arithmetic in hackage http://hackage.haskell.org/package/intervals-0.2.0 However it does not currently properly adjust the floating point rounding mode so containment isn't perfect. However, we are actively working on fixing up the Haskell MPFR bindings, which will let us reliably set rounding modes, cleaning up the interval arithmetic library to be just a little bit more pedantic. Due to the way GHC interacts with GMP this is a disturbingly difficult process. I have an unreleased library for working with Taylor models that builds on top of that and my automatic differentiation library, but without working MPFR bindings, it isn't sufficiently accurate for me to comfortably release. -Edward

On Sun, Mar 20, 2011 at 4:27 PM, Edward Amsden wrote:

...
Hi cafe,

I'm looking for a library that provides an instance of Num, Fractional, Floating, etc, but carries uncertainty values through calculations. A scan of hackage didn't turn anything up. Does anyone know of a library like this?

Thanks!

-- Edward Amsden Student Computer Science Rochester Institute of Technology www.edwardamsden.com

_______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe

_______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe

-- Edward Amsden Student Computer Science Rochester Institute of Technology www.edwardamsden.com

Tom Nielsen

10:47 p.m.

so if you want to do it the quick, easy and approximately correct, why not just use a monte-carlo monad? say you have two values x and y, for which you know the means and standard deviations mx, my, sdx, sdy. and you have a function f(x,y) that expresses wheat you want to know. so then in the random-fu monad you say (i can't member exactly how it works; the following is valid in "probably", my unreleased probability/stats library at https://github.com/glutamate/samfun) xs = sampleN 1000 $ do x <- normal mx sdx y <- normal my sdy return $ f(x,y) and then take the mean and standard deviation of xs? that's every bit as correct as propagating the uncertainty of f, except for the finite number of samples. (assuming your original uncertainties are gaussian) Tom On Sun, Mar 20, 2011 at 11:59 PM, Edward Amsden wrote:

...

I'm actually a CS undergrad in a physics lab class. I have permission from my professor to use computer programs for analysis of lab data. I need to do calculations on data with uncertainty, but uncertainty analysis on many formulae in physics is rather tedious. I was hoping for something with instances for Num, Fractional, Floating, etc. that would allow me to combine two uncertain values and get a new value with uncertainty. I've been working on writing one myself and I don't find the concept hard, but it's a lot of effort that I don't want to duplicate if it's been done already.

On Sun, Mar 20, 2011 at 5:46 PM, Tom Nielsen wrote:

...
Interval arithmetic is of course not the same as uncertainty, although computer scientists like to pretend that is the case. (and uncertainty estimates do not have the be "rough".)

In general the propagation of errors depends on whether the errors are independent or not. The rules are given in Taylor: An introduction to Error analysis (1997). Interval artihmetic corresponds to the worst case of non-independent and non-random errors. In the case of independent of random errors, you get:

data Approximately a = a :+/-: a

instance Num a => Num (Approximately a) where (m1 :+/-: err1) + (m2 :+/-: err2) = (m1+m2) :+/-: (sqrt(err1^2+err2^2) (m1 :+/-: err1) - (m2 :+/-: err2) = (m1-m2) :+/-: (sqrt(err1^2+err2^2) (m1 :+/-: err1) * (m2 :+/-: err2) = (m1*m2) :+/-: (sqrt((err1/m1)^2+(err2/m2)^2)

the general rule is

if y = f xs where xs :: [Approximately a], i.e f :: [Approximately a] -> Approximately a

the error term= sqrt $ sum $ map (^2) $ map (\(ym :+/-: yerr) -> partial-derivative-of-yerr-with-respect-to-partial-ym * yerr) xs

You can verify these things by running your calculation through soem sort of randomness monad (monte-carlo or random-fu packages) Anyways, I ended up not going down this route this because probabilistic data analysis gives you the correct error estimate without propagating error terms.

Tom

PS if you're a scientist and your accuracy estimate is on the same order as your rounding error, your are doing pretty well :-) At least in my field...

On Sun, Mar 20, 2011 at 8:38 PM, Edward Kmett wrote:

...
I have a package for interval arithmetic in hackage http://hackage.haskell.org/package/intervals-0.2.0 However it does not currently properly adjust the floating point rounding mode so containment isn't perfect. However, we are actively working on fixing up the Haskell MPFR bindings, which will let us reliably set rounding modes, cleaning up the interval arithmetic library to be just a little bit more pedantic. Due to the way GHC interacts with GMP this is a disturbingly difficult process. I have an unreleased library for working with Taylor models that builds on top of that and my automatic differentiation library, but without working MPFR bindings, it isn't sufficiently accurate for me to comfortably release. -Edward

On Sun, Mar 20, 2011 at 4:27 PM, Edward Amsden wrote:

...
Hi cafe,

I'm looking for a library that provides an instance of Num, Fractional, Floating, etc, but carries uncertainty values through calculations. A scan of hackage didn't turn anything up. Does anyone know of a library like this?

Thanks!

-- Edward Amsden Student Computer Science Rochester Institute of Technology www.edwardamsden.com

_______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe

_______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe

-- Edward Amsden Student Computer Science Rochester Institute of Technology www.edwardamsden.com

_______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe

Edward Kmett

21 Mar 21 Mar

1:36 p.m.

On Sun, Mar 20, 2011 at 5:46 PM, Tom Nielsen wrote:

...

Interval arithmetic is of course not the same as uncertainty, although computer scientists like to pretend that is the case. (and uncertainty estimates do not have the be "rough".)

Very true.

...

In general the propagation of errors depends on whether the errors are independent or not. The rules are given in Taylor: An introduction to Error analysis (1997). Interval artihmetic corresponds to the worst case of non-independent and non-random errors. In the case of independent of random errors, you get:

...

data Approximately a = a :+/-: a

instance Num a => Num (Approximately a) where (m1 :+/-: err1) + (m2 :+/-: err2) = (m1+m2) :+/-: (sqrt(err1^2+err2^2) (m1 :+/-: err1) - (m2 :+/-: err2) = (m1-m2) :+/-: (sqrt(err1^2+err2^2) (m1 :+/-: err1) * (m2 :+/-: err2) = (m1*m2) :+/-: (sqrt((err1/m1)^2+(err2/m2)^2)

the general rule is

if y = f xs where xs :: [Approximately a], i.e f :: [Approximately a] -> Approximately a

the error term= sqrt $ sum $ map (^2) $ map (\(ym :+/-: yerr) -> partial-derivative-of-yerr-with-respect-to-partial-ym * yerr) xs

The danger here is, of course, the side-condition of independence, which can make inhabitants of that type very difficult to reason about. e.g. x + x and 2*x in that world are very different. In this sense the interval arithmetic bounds _are_ safer to work with in the absence of sharing information even if they are less useful in the hands of an expert.

...

You can verify these things by running your calculation through soem sort of randomness monad (monte-carlo or random-fu packages) Anyways, I ended up not going down this route this because probabilistic data analysis gives you the correct error estimate without propagating error terms.

We are in total agreement here. =)

...

Tom

PS if you're a scientist and your accuracy estimate is on the same order as your rounding error, your are doing pretty well :-) At least in my field...

True enough, but in the case of interval arithmetic I like to be able to preserve the invariant that if I am working with intervals (even if only to collect accumulated rounding error in a Taylor model) that the answer lies within the interval, and doesn't escape due to some tight boundary condition or accumulated rounding error from when I was working too close to a pole. In the case of Taylor models we try to keep the size of the intervals as small as possible by using the first k terms of a Taylor polynomial and only catching the slop in an interval.This is important because of course adding and multiplying intervals will cause the size of the intervals to baloon very quickly. Since the intervals in question are very close to the scale of floating point rounding error as possible, and we often have to conservatively slop the rounding error over from the Taylor coefficients into the interval, accurate handling of tight corner cases is critical. -Edward

Edward Amsden

2:42 p.m.

...

...
PS if you're a scientist and your accuracy estimate is on the same order as your rounding error, your are doing pretty well :-) At least in my field...

True enough, but in the case of interval arithmetic I like to be able to preserve the invariant that if I am working with intervals (even if only to collect accumulated rounding error in a Taylor model) that the answer lies within the interval, and doesn't escape due to some tight boundary condition or accumulated rounding error from when I was working too close to a pole.

In the case of Taylor models we try to keep the size of the intervals as small as possible by using the first k terms of a Taylor polynomial and only catching the slop in an interval.This is important because of course adding and multiplying intervals will cause the size of the intervals to baloon very quickly. Since the intervals in question are very close to the scale of floating point rounding error as possible, and we often have to conservatively slop the rounding error over from the Taylor coefficients into the interval, accurate handling of tight corner cases is critical.

-Edward

So I'm feeling a bit elated that I've sparked my first theoretical discussion in cafe, though I don't have much to contribute. :\ However in the interests of the original question, I guess I should clarify. What we do in our physics class seems to be what is being called "interval analysis" in this discussion. We have experimental values with absolute uncertainties, and we need to propagate those uncertainties in a deterministic way through formulas. I don't think my professor would take kindly to a random sampling approach. The intervals library seemed a bit like what I'm looking for, except that it appears to be broken for the later ghc 6 versions and ghc 7. -- Edward Amsden Student Computer Science Rochester Institute of Technology www.edwardamsden.com

Edward Kmett

4:57 p.m.

On Mon, Mar 21, 2011 at 2:42 PM, Edward Amsden wrote:

...

So I'm feeling a bit elated that I've sparked my first theoretical discussion in cafe, though I don't have much to contribute. :\

However in the interests of the original question, I guess I should clarify.

What we do in our physics class seems to be what is being called "interval analysis" in this discussion. We have experimental values with absolute uncertainties, and we need to propagate those uncertainties in a deterministic way through formulas. I don't think my professor would take kindly to a random sampling approach.

The intervals library seemed a bit like what I'm looking for, except that it appears to be broken for the later ghc 6 versions and ghc 7.

The package should build fine, but hackage was flipping out because I commented out a pattern guard, and it looked like a misplaced haddock comment. I've pushed a new version of intervals to mollify hackage. It (or the old version) should cabal install just fine. -Edward Kmett

Carter Schonwald

6:01 p.m.

by the way, the link to the patch-tag repo for your intervals lib seems to be dead / patch-tag gets confused, is it that the link is outdated or that there are problems on patch-tag? -Carter On Mon, Mar 21, 2011 at 4:57 PM, Edward Kmett wrote:

...

On Mon, Mar 21, 2011 at 2:42 PM, Edward Amsden wrote:

...
So I'm feeling a bit elated that I've sparked my first theoretical discussion in cafe, though I don't have much to contribute. :\

However in the interests of the original question, I guess I should clarify.

What we do in our physics class seems to be what is being called "interval analysis" in this discussion. We have experimental values with absolute uncertainties, and we need to propagate those uncertainties in a deterministic way through formulas. I don't think my professor would take kindly to a random sampling approach.

The intervals library seemed a bit like what I'm looking for, except that it appears to be broken for the later ghc 6 versions and ghc 7.

The package should build fine, but hackage was flipping out because I commented out a pattern guard, and it looked like a misplaced haddock comment.

I've pushed a new version of intervals to mollify hackage.

It (or the old version) should cabal install just fine.

-Edward Kmett

_______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe

Carter Schonwald

6:02 p.m.

i'm now seeing that they've been moved to github, never mind! On Mon, Mar 21, 2011 at 6:01 PM, Carter Schonwald < carter.schonwald@gmail.com> wrote:

...

by the way, the link to the patch-tag repo for your intervals lib seems to be dead / patch-tag gets confused, is it that the link is outdated or that there are problems on patch-tag?

-Carter

On Mon, Mar 21, 2011 at 4:57 PM, Edward Kmett wrote:

...
On Mon, Mar 21, 2011 at 2:42 PM, Edward Amsden wrote:

...
So I'm feeling a bit elated that I've sparked my first theoretical discussion in cafe, though I don't have much to contribute. :\

However in the interests of the original question, I guess I should clarify.

What we do in our physics class seems to be what is being called "interval analysis" in this discussion. We have experimental values with absolute uncertainties, and we need to propagate those uncertainties in a deterministic way through formulas. I don't think my professor would take kindly to a random sampling approach.

The intervals library seemed a bit like what I'm looking for, except that it appears to be broken for the later ghc 6 versions and ghc 7.

The package should build fine, but hackage was flipping out because I commented out a pattern guard, and it looked like a misplaced haddock comment.

I've pushed a new version of intervals to mollify hackage.

It (or the old version) should cabal install just fine.

-Edward Kmett

_______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe

Edward Kmett

22 Mar 22 Mar

3:30 p.m.

I moved all of my repositories over to github back around June: http://twitter.com/#!/kmett/status/16174477854 I should update that link. -Edward On Mon, Mar 21, 2011 at 6:02 PM, Carter Schonwald < carter.schonwald@gmail.com> wrote:

...

i'm now seeing that they've been moved to github, never mind!

On Mon, Mar 21, 2011 at 6:01 PM, Carter Schonwald < carter.schonwald@gmail.com> wrote:

...
by the way, the link to the patch-tag repo for your intervals lib seems to be dead / patch-tag gets confused, is it that the link is outdated or that there are problems on patch-tag?

-Carter

On Mon, Mar 21, 2011 at 4:57 PM, Edward Kmett wrote:

...
On Mon, Mar 21, 2011 at 2:42 PM, Edward Amsden wrote:

...
So I'm feeling a bit elated that I've sparked my first theoretical discussion in cafe, though I don't have much to contribute. :\

However in the interests of the original question, I guess I should clarify.

What we do in our physics class seems to be what is being called "interval analysis" in this discussion. We have experimental values with absolute uncertainties, and we need to propagate those uncertainties in a deterministic way through formulas. I don't think my professor would take kindly to a random sampling approach.

The intervals library seemed a bit like what I'm looking for, except that it appears to be broken for the later ghc 6 versions and ghc 7.

The package should build fine, but hackage was flipping out because I commented out a pattern guard, and it looked like a misplaced haddock comment.

I've pushed a new version of intervals to mollify hackage.

It (or the old version) should cabal install just fine.

-Edward Kmett

_______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe

Tom Nielsen

21 Mar 21 Mar

5:41 p.m.

...

The danger here is, of course, the side-condition of independence, which can make inhabitants of that type very difficult to reason about. e.g. x + x and 2*x in that world are very different.

Yes. I was surprised (maybe i shouldn't have been): sampler = do x <- gauss 0.0 1.0 y <- gauss 0.0 1.0 return $ (2*x, x+y) main = do xys <- take 100000 `fmap` runSamplerIO sampler print $ runStat (both (before varF fst) (before varF snd)) $ xys => (3.9988971177326498,2.0112123117664975) that is, (variance of 2*x, variance of x+y) The problem is that when you say x+x you don't really mean it; you mean something like liftM2(+) xdist xdist in a probability monad. Had I changed "x+y" to "x+x", I would obviously have gotten identical variances. So maybe referential transparency is not lost after all. Tom

Henning Thielemann

6:24 p.m.

On Mon, 21 Mar 2011, Tom Nielsen wrote:

...

sampler = do x <- gauss 0.0 1.0 y <- gauss 0.0 1.0 return $ (2*x, x+y)

main = do xys <- take 100000 `fmap` runSamplerIO sampler print $ runStat (both (before varF fst) (before varF snd)) $ xys

=> (3.9988971177326498,2.0112123117664975) that is, (variance of 2*x, variance of x+y)

Variances of independent random variables (here x and y) are additive, for dependent variables (say, x and x) this does not hold.

5175

Age (days ago)

5177

Last active (days ago)

List overview

Download

13 comments

5 participants

participants (5)

Carter Schonwald
Edward Amsden
Edward Kmett
Henning Thielemann
Tom Nielsen