0

I have to calculate the probability of the proper decoding of a bit copied n times. The following formula is supposed to be the answer:

the formula

In Haskell, I coded it as follows:

fac 1 = 1
fac n = fac (n-1) * n
--prob :: (Integral b, Fractional b) => (b, b) -> b
--prob :: (Int,Int) -> Double
prob (n, k)
    | n==k = (0.01**k)
    | otherwise = factor (n, k) * (0.01 ** k) * (0.99**(n-k)) + prob (n, (k+1)) 
    where
        factor (n, k) = (fac n / ((fac k)* (fac n-k)))

1 - prob (3,2) gives the result 0.99992575, which is incorrect, as it should be 0.99970. Does anyone know where I went wrong?

Mark Karpov
  • 7,499
  • 2
  • 27
  • 62
Dirk Huethorst
  • 1
  • 1
    Make the epsilon a parameter to `prob`, and I personally would make a wrapper around it that calculates `k` for me from the `n` passed in. Adjust the epsilon to be smaller and see if that's what the problem is. – bheklilr Sep 25 '14 at 14:29
  • why would you comment out the type? – user3329719 Sep 25 '14 at 15:14
  • this looks too much like newtons binom formula – user3329719 Sep 25 '14 at 15:15
  • @user3329719 They're very related: This formula calculates the probability of `(n+1)/2` or more errors using a [binomial distribution](https://en.wikipedia.org/wiki/Binomial_distribution). – Ørjan Johansen Sep 26 '14 at 06:32

2 Answers2

6

The reason is function precedence. if you would look into the definition of prob you will see this:

(fac n-k)

because function application has the most precedence this is parsed as

((fac n) - k)

so your code should be

(fac (n-k))

which gives a result of 0.999702 on my computer.

user3329719
  • 230
  • 1
  • 8
0

These are a couple of best-practices the code is lacking. I infact already answered the question itself.

1- do not use tuples as input. in Haskell, functions can have multiple parameters. the syntax is f x y for calling f on x and y. the types also have similar syntax. this transforms your code to:

fac 1 = 1
fac n = fac (n-1) * n
--prob :: (Integral b, Fractional b) => b -> b -> b (two parameters of type b and output of type b)
--prob :: Int -> Int -> Double
prob n k
    | n==k = (0.01**k)
    | otherwise = factor n k * (0.01 ** k) * (0.99**(n-k)) + prob n (k+1) 
    where
        factor n k = (fac n / ((fac k)* (fac (n-k))))

2- if you will notice, fac will only work on Integers, and similarly does factor. prob infact then has type of (Fractional a, Integral b) -> b -> b -> a or alternatively Integer -> Integer -> Float. why not give them their true type?

this transformation requires changing ** (which gets two floating point numbers) to ^ (which gets an Integer as it's second parameter) and using the function fromIntegral which casts an Integer to an arbitrary number Data. 

fac :: Integral a => a -> a -- or alternatively Integer -> Integer
fac 1 = 1
fac n = fac (n-1) * n
prob n k
    | n==k = (0.01 ^^ k)
    | otherwise = fromIntegral (factor n k) * (0.01 ^^ k) * (0.99 ^^ (n-k) + prob n (k+1) 
    where
        factor n k = div (fac n) (fac k * fac (n-k)) -- div is Integer division operator.

now prob has type of (Integral a, Floating b) => a -> a -> b which means it gets two parameters of type a (which is an Integral instance) and returns a value of type b.

user3329719
  • 230
  • 1
  • 8