使用 Haskell 'statistics' 包进行线性回归

Question

我需要对一组测量数据进行线性回归。我知道statistics 是去图书馆工作。 Statistics.Regression 模块中ols 函数的描述表明它是我需要的东西。但是，从文档中不清楚如何表示我的测量数据，例如

measurements = [(1.0, 2.0), (2.0, 2.5), (3.0, 3.0)] :: [(Double, Double)]  -- (X, Y) coordinates of points

到签名为的函数

ols :: Matrix   -- A has at least as many rows as columns.
    -> Vector   -- b has the same length as columns in A.
    -> Vector

以及如何解释结果。 我希望得到等式y = a + bx 的值a = 1.5 和b = 0.5。

我必须如何构造输入 Matrix 和 Vector，以及生成的 Vector 中的内容是什么？

Answer 1

以某种方式期望统计库提供一些通用性的功能，因此如果您追求一个简单的案例，您必须缩小提供的功能。在这里，库函数需要几个因果变量，而不仅仅是一个。此外，它使用向量和矩阵，您可能只需要简单的列表。

为了展示如何使用带有单个因果变量的库，我假设您想要一个只使用普通 Haskell 数据类型的接口，如下所示：

-- [(x0, y0), (x1, y1), , (x2, y2) ... ] -> ((a, b), r2)  for  y ≃ a + b*x
simpleRegression1 :: [(Double, Double)] -> ((Double, Double), Double)

使用olsRegress 似乎比使用ols 更简单，因为：

1. 不需要矩阵数据类型
2. 您可以轻松获得 y 截距值（即表示法中的 a）
3. 您还可以从同一个调用中获得拟合优度系数

代码可以写成如下：

import qualified  Data.List              as  L
import qualified  Data.Vector.Unboxed    as  DVU
import qualified  Statistics.Regression  as  SR

-- [(x0, y0), (x1, y1), , (x2, y2) ... ] -> ((a, b), r2)  for  y ≃ a + b*x
simpleRegression1 :: [(Double, Double)] -> ((Double, Double), Double)
simpleRegression1 xyPairs =
   let  xList    = L.map  fst  xyPairs
        yList    = L.map  snd  xyPairs
        xVecList = [DVU.fromList  xList]
        yVec     =  DVU.fromList  yList
        (sv, r2) = SR.olsRegress xVecList yVec
        [b, a]   = DVU.toList sv
    in
        ((a,b), r2)

在我们的例子中，只有一个列向量，所以xVecList 只有一个元素。 正如olsRegress documentation 中提到的，y 截距“a”值是输出向量的最后一个元素。

测试代码：

main = do
    let  measurements = [(1.0, 1.9999), (2.0, 2.5001), (3.0, 2.9999)]
         ---- measurements = [(1.0, 2.0), (2.0, 2.5), (3.0, 3.0)]
         ((a,b), r2)  = simpleRegression1  measurements
    putStrLn $ "measurements = " ++ (show measurements)
    putStrLn $ "a = " ++ (show a) ++ "  b = " ++ (show b) ++
               "  r2 = " ++ (show r2)

    let  ypList   = L.map  (\x -> a+b*x)  (L.map fst measurements)
         diffList = L.zipWith  (-)  (L.map snd measurements)  ypList
    putStrLn $ "ypList = " ++ (show ypList)
    putStrLn $ "diffList = " ++ (show diffList)

测试输出：

measurements = [(1.0,1.9999),(2.0,2.5001),(3.0,2.9999)]
a = 1.4999666666666656  b = 0.5000000000000006  r2 = 0.9999999466666695
ypList = [1.999966666666666,2.4999666666666664,2.9999666666666673]
diffList = [-6.666666666599319e-5,1.3333333333376274e-4,-6.66666666675475e-5]

注意还有一个LinearRegression module。

Answer 2

让我重写ols 的签名以反映它在概念上实际完成的工作：

ols :: (VectorSpace u, VectorSpace v) => (u +> v) -> v -> u

其中+> 表示linear function，这就是矩阵所代表的内容。

所以，u 是您想要的结果（即参数a 和b）。矩阵参数是将这些参数映射到实际测量的预测的（线性）函数。因此，这涉及在数据集中的每个 x 位置计算 a·x + b。

在更受数学启发（而不是受 Matlab 或 R 启发）的库中，您可以大致这样编写它：

m :: Matrix
m = fromFunction $ \(a,b) -> fmap (\x -> a*x + b) [1,2,3])

然后v 类型的测量向量是您在输出中拥有的y 值，因此您将评估

ols m [2, 2.5, 3]

事实上，在linearmap-category 中——与statistics 不同的是，它的类型是正确的——你几乎可以通过这种方式做到这一点：

> :m +@987654323@
> import @987654324@ (V2(..), V3(..))
> lfun (\(V2 a b) -> fmap (\x -> a*x + b) (V3 1 2 3)) @987654325@ V3 2 2.5 3
V2 0.4999999999999998 1.5000000000000018

在statisics 中，我怀疑您自己需要从线性函数中构建矩阵。这有点繁琐但并不困难：构建矩阵意味着简单地放入可能的基础输入(1,0) 和(0,1)，为它们评估函数（即一次用于a=1，b=0，一次用于a=0， b=1) 并将所有x 的结果记录为列向量。