什么相当于Haskell中Ruby的pnormaldist统计函数？

Question

如下所示:http://www.evanmiller.org/how-not-to-sort-by-average-rating.html

这是Ruby代码本身,在Statistics2库中实现:

# inverse of normal distribution ([2])
# Pr( (-\infty, x] ) = qn -> x
def pnormaldist(qn)
  b = [1.570796288, 0.03706987906, -0.8364353589e-3,
       -0.2250947176e-3, 0.6841218299e-5, 0.5824238515e-5,
       -0.104527497e-5, 0.8360937017e-7, -0.3231081277e-8,
       0.3657763036e-10, 0.6936233982e-12]

  if(qn < 0.0 || 1.0 < qn)
    $stderr.printf("Error : qn <= 0 or qn >= 1  in pnorm()!\n")
    return 0.0;
  end
  qn == 0.5 and return 0.0

  w1 = qn
  qn > 0.5 and w1 = 1.0 - w1
  w3 = -Math.log(4.0 * w1 * (1.0 - w1))
  w1 = b[0]
  1.upto 10 do |i|
    w1 += b[i] * w3**i;
  end
  qn > 0.5 and return Math.sqrt(w1 * w3)
  -Math.sqrt(w1 * w3)
end

Answer 1

在Hackage上挖掘,有许多统计库:

• hmatrix-gsl-stats - 与GSL的纯粹绑定
• hstatistics - GSL的更高级别接口
• hstats - 常见的统计方法
• 统计 - 更常见的统计方法
• statistics-linreg - 两个样本之间的线性回归,基于其他统计包.

你想要一个版本pnormaldist,"返回normaldist(x)的P值".

• Statistics.Distribution.Normal,来自统计包,提供了许多用于操作正态分布的函数.
• Statistics.Test.NonParametric包含许多与P值有关的事情.

也许某些东西提供了你需要的东西？

Answer 2

翻译非常简单:

module PNormalDist where

pnormaldist :: (Ord a, Floating a) => a -> Either String a
pnormaldist qn
  | qn < 0 || 1 < qn = Left "Error: qn must be in [0,1]"
  | qn == 0.5        = Right 0.0
  | otherwise        = Right $
      let w3 = negate . log $ 4 * qn * (1 - qn)
          b = [ 1.570796288, 0.03706987906, -0.8364353589e-3, 
                -0.2250947176e-3, 0.6841218299e-5, 0.5824238515e-5, 
                -0.104527497e-5, 0.8360937017e-7, -0.3231081277e-8, 
                0.3657763036e-10, 0.6936233982e-12]
          w1 = sum . zipWith (*) b $ iterate (*w3) 1
      in (signum $ qn - 0.5) * sqrt (w1 * w3)

首先,让我们看一下ruby - 它会返回一个值,但有时它会输出一条错误消息(当给出一个不正确的参数时).这不是很糟糕,所以让我们的返回值为Either String a- Left String如果给出不正确的参数,我们将返回带错误消息的地方,Right a否则.

现在我们检查顶部的两个案例:

• qn < 0 || 1 < qn = Left "Error: qn must be in [0,1]"- 这是错误条件,qn超出范围时.
• qn == 0.5 = Right 0.0 - 这是红宝石检查 qn == 0.5 and return * 0.0

接下来,我们w1在ruby代码中定义.但是我们稍后重新定义了几行,这不是非常红宝石.我们w1第一次存储的值会立即用于定义中w3,为什么我们不跳过存储它w1？我们甚至不需要这样做qn > 0.5 and w1 = 1.0 - w1,因为我们w1 * (1.0 - w1)在w3的定义中使用了该产品.

所以我们跳过所有这些,直接进入定义w3 = negate . log $ 4 * qn * (1 - qn).

接下来是定义b,它是ruby代码的直接提升(ruby的数组文字语法是haskell的列表语法).

这是最棘手的一点 - 定义最终价值w3.ruby代码的作用

w1 = b[0]
1.upto 10 do |i|
  w1 += b[i] * w3**i;
end

是什么称为折叠 - 将一组值(存储在ruby数组中)减少为单个值.我们可以使用以下函数更多地重述(但仍然在ruby中)Array#reduce:

w1 = b.zip(0..10).reduce(0) do |accum, (bval,i)|
  accum + bval * w3^i
end

注意我是如何b[0]使用身份进入循环的b[0] == b[0] * w3^0.

现在我们可以直接将它移植到haskell,但它有点难看

w1 = foldl 0 (\accum (bval,i) -> accum + bval * w3**i) $ zip b [0..10]

相反,我把它分成了几个步骤 - 首先,我们并不真正需要i,我们只需要w3(从开始w3^0 == 1)的权力,所以让我们计算那些iterate (*w3) 1.

然后,我们最终只需要他们的产品,而不是使用b的元素将它们拉成一对,因此我们可以使用它们将它们压缩到每对产品中zipWith (*) b.

现在我们的折叠功能非常简单 - 我们只需要总结我们可以使用的产品sum.

最后,sqrt (w1 * w3)根据qn是大于还是小于0.5(我们已经知道它不相等),我们决定是否返回正负.因此,我不是像在ruby代码中那样在两个不同的位置计算平方根,而是计算它一次,并将其乘以+1或-1根据符号qn - 0.5(signum 只返回值的符号).