statistics-0.16.5.0: tests/Tests/Parametric.hs
module Tests.Parametric (tests) where
import Data.Maybe (fromJust)
import Statistics.Test.StudentT
import Statistics.Types
import qualified Data.Vector.Unboxed as U
import qualified Data.Vector as V
import Test.Tasty (testGroup, TestTree)
import Test.Tasty.HUnit (testCase, assertBool)
import Tests.Helpers (testEquality)
import qualified Test.Tasty as Tst
import Statistics.Test.Levene
import Statistics.Test.Bartlett
tests :: Tst.TestTree
tests = testGroup "Parametric tests" [studentTTests, bartlettTests, leveneTests]
-- 2 samples x 20 obs data
--
-- Both samples are samples from normal distributions with the same variance (= 1.0),
-- but their means are different (0.0 and 0.5, respectively).
--
-- You can reproduce the data with R (3.1.0) as follows:
-- set.seed(0)
-- sample1 = rnorm(20)
-- sample2 = rnorm(20, 0.5)
-- student = t.test(sample1, sample2, var.equal=T)
-- welch = t.test(sample1, sample2)
-- paired = t.test(sample1, sample2, paired=T)
sample1, sample2 :: U.Vector Double
sample1 = U.fromList [
1.262954284880793e+00,
-3.262333607056494e-01,
1.329799262922501e+00,
1.272429321429405e+00,
4.146414344564082e-01,
-1.539950041903710e+00,
-9.285670347135381e-01,
-2.947204467905602e-01,
-5.767172747536955e-03,
2.404653388857951e+00,
7.635934611404596e-01,
-7.990092489893682e-01,
-1.147657009236351e+00,
-2.894615736882233e-01,
-2.992151178973161e-01,
-4.115108327950670e-01,
2.522234481561323e-01,
-8.919211272845686e-01,
4.356832993557186e-01,
-1.237538421929958e+00]
sample2 = U.fromList [
2.757321147216907e-01,
8.773956459817011e-01,
6.333363608148415e-01,
1.304189509744908e+00,
4.428932256161913e-01,
1.003607972233726e+00,
1.585769362145687e+00,
-1.909538396968303e-01,
-7.845993538721883e-01,
5.467261721883520e-01,
2.642934435604988e-01,
-4.288825501025439e-02,
6.668968254321778e-02,
-1.494716467962331e-01,
1.226750747385451e+00,
1.651911754087200e+00,
1.492160365445798e+00,
7.048689050811874e-02,
1.738304100853380e+00,
2.206537181457307e-01]
testTTest :: String
-> PValue Double
-> Test d
-> [Tst.TestTree]
testTTest name pVal test =
[ testEquality name (isSignificant pVal test) NotSignificant
, testEquality name (isSignificant (mkPValue $ pValue pVal + 1e-5) test)
Significant
]
studentTTests :: Tst.TestTree
studentTTests = testGroup "StudentT test" $ concat
[ -- R: t.test(sample1, sample2, alt="two.sided", var.equal=T)
testTTest "two-sample t-test SamplesDiffer Student"
(mkPValue 0.03410) (fromJust $ studentTTest SamplesDiffer sample1 sample2)
-- R: t.test(sample1, sample2, alt="two.sided", var.equal=F)
, testTTest "two-sample t-test SamplesDiffer Welch"
(mkPValue 0.03483) (fromJust $ welchTTest SamplesDiffer sample1 sample2)
-- R: t.test(sample1, sample2, alt="two.sided", paired=T)
, testTTest "two-sample t-test SamplesDiffer Paired"
(mkPValue 0.03411) (fromJust $ pairedTTest SamplesDiffer sample12)
-- R: t.test(sample1, sample2, alt="less", var.equal=T)
, testTTest "two-sample t-test BGreater Student"
(mkPValue 0.01705) (fromJust $ studentTTest BGreater sample1 sample2)
-- R: t.test(sample1, sample2, alt="less", var.equal=F)
, testTTest "two-sample t-test BGreater Welch"
(mkPValue 0.01741) (fromJust $ welchTTest BGreater sample1 sample2)
-- R: t.test(sample1, sample2, alt="less", paired=F)
, testTTest "two-sample t-test BGreater Paired"
(mkPValue 0.01705) (fromJust $ pairedTTest BGreater sample12)
]
where sample12 = U.zip sample1 sample2
------------------------------------------------------------
-- Bartlett's Test
------------------------------------------------------------
bartlettTests :: TestTree
bartlettTests = testGroup "Bartlett's test"
[ testCase "a,b,c" $ testBartlettTest [a,b,c] 1.8027132567760222 0.40601846976301237
, testCase "a,b" $ testBartlettTest [a,b] 0.005221063776321886 0.9423974408021293
, testCase "a,c" $ testBartlettTest [a,c] 1.1531619271845452 0.2828882244527482
, testCase "a,a" $ testBartlettTest [a,a] 0.0 1.0
]
where
a = U.fromList [9.88, 9.12, 9.04, 8.98, 9.00, 9.08, 9.01, 8.85, 9.06, 8.99]
b = U.fromList [8.88, 8.95, 9.29, 9.44, 9.15, 9.58, 9.36, 9.18, 8.67, 9.05]
c = U.fromList [8.95, 8.12, 8.95, 8.85, 8.03, 8.84, 8.07, 8.98, 8.86, 8.98]
testBartlettTest
:: [U.Vector Double]
-> Double
-> Double
-> IO ()
testBartlettTest samples w p = do
r <- case bartlettTest samples of
Left _ -> error "Bartlett's test failed"
Right r -> pure r
approxEqual "W" 1e-9 (testStatistics r) w
approxEqual "p" 1e-9 (pValue $ testSignificance r) p
------------------------------------------------------------
-- Levene's Test (Trimmed Mean)
------------------------------------------------------------
leveneTests :: TestTree
leveneTests = testGroup "Levene test"
-- Statistics' value and p-values are computed using
[ testCase "a,b,c Mean" $ testLeveneTest [a,b,c] Mean 7.905194483442054 0.001983795817472731
, testCase "a,b Mean" $ testLeveneTest [a,b] Mean 8.83873787256358 0.008149720958328811
, testCase "a,a Mean" $ testLeveneTest [a,a] Mean 0.0 1.0
, testCase "a,b,c Median" $ testLeveneTest [a,b,c] Median 7.584952754501659 0.002431505967249681
, testCase "a,b Median" $ testLeveneTest [a,b] Median 8.461374333228711 0.009364737715584399
, testCase "aL,bL Mean" $ testLeveneTest [aL,bL] Mean 5.84424549939465 0.01653410652558999
, testCase "aL,bL Trimmed" $ testLeveneTest [aL,bL] (Trimmed 0.05) 8.368311226366314 0.004294953946529551
]
where
a = V.fromList [8.88, 9.12, 9.04, 8.98, 9.00, 9.08, 9.01, 8.85, 9.06, 8.99]
b = V.fromList [8.88, 8.95, 9.29, 9.44, 9.15, 9.58, 8.36, 9.18, 8.67, 9.05]
c = V.fromList [8.95, 9.12, 8.95, 8.85, 9.03, 8.84, 9.07, 8.98, 8.86, 8.98]
-- Large samples for testing trimmed
aL = V.fromList [
-0.18919252, -1.62837673, 5.21332355, -0.00962043, -0.28417847,
-0.88128233, 1.49698436, 6.1780359 , -1.22301348, 3.34598245,
5.33227264, -0.88732069, 0.14487346, 2.61060215, 4.22033907,
2.53139215, -0.72131061, 0.53063607, -0.60510374, -0.73230842,
1.54037043, -2.81103963, 3.40763063, 0.49005324, 2.13085513,
5.68650547, 4.16397279, -0.17325097, 1.12664972, 4.23297516,
4.15943436, -1.01452078, 2.40391646, 0.83019962, 0.29665879,
-3.83031046, -1.98576933, 1.5356527 , 1.30773365, 0.292818 ,
2.45877828, 1.06482289, -0.63241873, 1.58465379, 1.96577614,
2.25791943, 4.13769848, -2.38595767, -0.65801423, -2.54007791,
3.17428087, 4.32096964, 0.92240335, -2.38101319, 1.35692587,
1.48279101, -0.04438309, 0.50296642, 2.08261495, 1.33181215,
-1.95427198, 4.95406809, 1.51294898, -2.68536129, -0.2441218 ,
2.41142613, 4.71051493, 2.66618697, 1.12668301, -0.25732583,
1.25021838, -1.27523641, 5.01638744, 3.38864442, 0.17979744,
-0.88481645, 3.89346357, -0.51512217, -1.60542888, 0.88378679,
-2.12962732, -1.35989539, 5.09215112, -1.37442481, 0.83578405,
0.13829571, 1.25171481, 3.60552158, -3.24051591, -0.44301834,
0.78253445, 1.76098254, 1.79677434, -0.19010505, 3.07640466,
3.02853882, 1.24849063, 4.84505382, 6.82274999, 2.24063474]
bL = V.fromList [
2.15584101, -2.74876744, -0.82231894, 1.97518087, 2.59280595,
1.28703417, 2.40450278, 1.9761031 , 2.35186598, 1.15611047,
2.26709318, 1.2832138 , -2.1486074 , 0.27563011, -0.51816861,
0.89658424, 3.27069545, 1.72846646, 3.84454277, 5.58301459,
-0.40878188, 3.41602853, 1.1281526 , 0.9665913 , 0.76567084,
1.69522855, 1.69133014, 0.70529264, 2.65243202, -1.0088019 ,
-0.62431026, 3.76667396, 3.66225181, 0.73217579, 0.04478736,
0.4169833 , 0.77065631, -1.31484093, 1.23858618, -0.08339456,
3.14154286, 1.84358218, -0.53511423, -3.4919477 , 0.24076997,
3.59381684, 1.99497806, 2.95499775, 1.67157731, 0.0214764 ,
3.32161612, -2.64762427, 0.06486472, 0.19653897, 1.34954235,
1.18568747, -0.54434597, -3.35544223, 1.41933109, 0.95100195,
2.7182116 , 1.1334068 , -0.95297806, -0.05421818, 1.42248799,
-3.96201277, -3.21309254, -0.21209211, 0.9689551 , 0.13526401,
-0.88656198, 0.41331783, -3.18766064, 4.34948246, 1.35656384,
0.41920101, -0.46578994, 1.55181583, 2.43937014, 2.49040644,
4.10505494, 1.68856296, 1.31503895, 0.41123368, 0.73242999,
0.2804349 , -1.83494592, -0.31073195, 2.61185513, 2.91645094,
1.26097638, 2.64197134, 3.88931972, 0.03783002, 2.55209729,
3.46869549, 0.96348003, 2.27658242, 2.7613171 , -0.1372434 ]
testLeveneTest
:: [V.Vector Double]
-> Center
-> Double
-> Double
-> IO ()
testLeveneTest samples center w p = do
r <- case levenesTest center samples of
Left _ -> error "Levene's test failed"
Right r -> pure r
approxEqual "W" 1e-9 (testStatistics r) w
approxEqual "p" 1e-9 (pValue $ testSignificance r) p
----------------------------------------------------------------
approxEqual :: String -> Double -> Double -> Double -> IO ()
approxEqual name epsilon actual expected =
assertBool (name ++ ": expected ≈ " ++ show expected ++ ", got " ++ show actual)
(diff < epsilon)
where
diff = abs (actual - expected)