{-# LANGUAGE OverloadedStrings #-}
{-# LANGUAGE DeriveDataTypeable, DeriveGeneric #-}
-- |
-- Module    : Statistics.Distribution.ChiSquared
-- Copyright : (c) 2010 Alexey Khudyakov
-- License   : BSD3
--
-- Maintainer  : bos@serpentine.com
-- Stability   : experimental
-- Portability : portable
--
-- The chi-squared distribution. This is a continuous probability
-- distribution of sum of squares of k independent standard normal
-- distributions. It's commonly used in statistical tests
module Statistics.Distribution.ChiSquared (
          ChiSquared
        , chiSquaredNDF
        -- * Constructors
        , chiSquared
        , chiSquaredE
        ) where

import Control.Applicative
import Data.Aeson            (FromJSON(..), ToJSON, Value(..), (.:))
import Data.Binary           (Binary(..))
import Data.Data             (Data, Typeable)
import GHC.Generics          (Generic)
import Numeric.SpecFunctions ( incompleteGamma,invIncompleteGamma,logGamma,digamma)
import Numeric.MathFunctions.Constants (m_neg_inf)
import qualified System.Random.MWC.Distributions as MWC

import qualified Statistics.Distribution         as D
import Statistics.Internal



-- | Chi-squared distribution
newtype ChiSquared = ChiSquared
  { ChiSquared -> Int
chiSquaredNDF :: Int
    -- ^ Get number of degrees of freedom
  }
  deriving (ChiSquared -> ChiSquared -> Bool
(ChiSquared -> ChiSquared -> Bool)
-> (ChiSquared -> ChiSquared -> Bool) -> Eq ChiSquared
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: ChiSquared -> ChiSquared -> Bool
$c/= :: ChiSquared -> ChiSquared -> Bool
== :: ChiSquared -> ChiSquared -> Bool
$c== :: ChiSquared -> ChiSquared -> Bool
Eq, Typeable, Typeable ChiSquared
DataType
Constr
Typeable ChiSquared
-> (forall (c :: * -> *).
    (forall d b. Data d => c (d -> b) -> d -> c b)
    -> (forall g. g -> c g) -> ChiSquared -> c ChiSquared)
-> (forall (c :: * -> *).
    (forall b r. Data b => c (b -> r) -> c r)
    -> (forall r. r -> c r) -> Constr -> c ChiSquared)
-> (ChiSquared -> Constr)
-> (ChiSquared -> DataType)
-> (forall (t :: * -> *) (c :: * -> *).
    Typeable t =>
    (forall d. Data d => c (t d)) -> Maybe (c ChiSquared))
-> (forall (t :: * -> * -> *) (c :: * -> *).
    Typeable t =>
    (forall d e. (Data d, Data e) => c (t d e))
    -> Maybe (c ChiSquared))
-> ((forall b. Data b => b -> b) -> ChiSquared -> ChiSquared)
-> (forall r r'.
    (r -> r' -> r)
    -> r -> (forall d. Data d => d -> r') -> ChiSquared -> r)
-> (forall r r'.
    (r' -> r -> r)
    -> r -> (forall d. Data d => d -> r') -> ChiSquared -> r)
-> (forall u. (forall d. Data d => d -> u) -> ChiSquared -> [u])
-> (forall u.
    Int -> (forall d. Data d => d -> u) -> ChiSquared -> u)
-> (forall (m :: * -> *).
    Monad m =>
    (forall d. Data d => d -> m d) -> ChiSquared -> m ChiSquared)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> ChiSquared -> m ChiSquared)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> ChiSquared -> m ChiSquared)
-> Data ChiSquared
ChiSquared -> DataType
ChiSquared -> Constr
(forall b. Data b => b -> b) -> ChiSquared -> ChiSquared
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> ChiSquared -> c ChiSquared
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c ChiSquared
forall a.
Typeable a
-> (forall (c :: * -> *).
    (forall d b. Data d => c (d -> b) -> d -> c b)
    -> (forall g. g -> c g) -> a -> c a)
-> (forall (c :: * -> *).
    (forall b r. Data b => c (b -> r) -> c r)
    -> (forall r. r -> c r) -> Constr -> c a)
-> (a -> Constr)
-> (a -> DataType)
-> (forall (t :: * -> *) (c :: * -> *).
    Typeable t =>
    (forall d. Data d => c (t d)) -> Maybe (c a))
-> (forall (t :: * -> * -> *) (c :: * -> *).
    Typeable t =>
    (forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c a))
-> ((forall b. Data b => b -> b) -> a -> a)
-> (forall r r'.
    (r -> r' -> r) -> r -> (forall d. Data d => d -> r') -> a -> r)
-> (forall r r'.
    (r' -> r -> r) -> r -> (forall d. Data d => d -> r') -> a -> r)
-> (forall u. (forall d. Data d => d -> u) -> a -> [u])
-> (forall u. Int -> (forall d. Data d => d -> u) -> a -> u)
-> (forall (m :: * -> *).
    Monad m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> Data a
forall u. Int -> (forall d. Data d => d -> u) -> ChiSquared -> u
forall u. (forall d. Data d => d -> u) -> ChiSquared -> [u]
forall r r'.
(r -> r' -> r)
-> r -> (forall d. Data d => d -> r') -> ChiSquared -> r
forall r r'.
(r' -> r -> r)
-> r -> (forall d. Data d => d -> r') -> ChiSquared -> r
forall (m :: * -> *).
Monad m =>
(forall d. Data d => d -> m d) -> ChiSquared -> m ChiSquared
forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d) -> ChiSquared -> m ChiSquared
forall (c :: * -> *).
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c ChiSquared
forall (c :: * -> *).
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> ChiSquared -> c ChiSquared
forall (t :: * -> *) (c :: * -> *).
Typeable t =>
(forall d. Data d => c (t d)) -> Maybe (c ChiSquared)
forall (t :: * -> * -> *) (c :: * -> *).
Typeable t =>
(forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c ChiSquared)
$cChiSquared :: Constr
$tChiSquared :: DataType
gmapMo :: (forall d. Data d => d -> m d) -> ChiSquared -> m ChiSquared
$cgmapMo :: forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d) -> ChiSquared -> m ChiSquared
gmapMp :: (forall d. Data d => d -> m d) -> ChiSquared -> m ChiSquared
$cgmapMp :: forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d) -> ChiSquared -> m ChiSquared
gmapM :: (forall d. Data d => d -> m d) -> ChiSquared -> m ChiSquared
$cgmapM :: forall (m :: * -> *).
Monad m =>
(forall d. Data d => d -> m d) -> ChiSquared -> m ChiSquared
gmapQi :: Int -> (forall d. Data d => d -> u) -> ChiSquared -> u
$cgmapQi :: forall u. Int -> (forall d. Data d => d -> u) -> ChiSquared -> u
gmapQ :: (forall d. Data d => d -> u) -> ChiSquared -> [u]
$cgmapQ :: forall u. (forall d. Data d => d -> u) -> ChiSquared -> [u]
gmapQr :: (r' -> r -> r)
-> r -> (forall d. Data d => d -> r') -> ChiSquared -> r
$cgmapQr :: forall r r'.
(r' -> r -> r)
-> r -> (forall d. Data d => d -> r') -> ChiSquared -> r
gmapQl :: (r -> r' -> r)
-> r -> (forall d. Data d => d -> r') -> ChiSquared -> r
$cgmapQl :: forall r r'.
(r -> r' -> r)
-> r -> (forall d. Data d => d -> r') -> ChiSquared -> r
gmapT :: (forall b. Data b => b -> b) -> ChiSquared -> ChiSquared
$cgmapT :: (forall b. Data b => b -> b) -> ChiSquared -> ChiSquared
dataCast2 :: (forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c ChiSquared)
$cdataCast2 :: forall (t :: * -> * -> *) (c :: * -> *).
Typeable t =>
(forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c ChiSquared)
dataCast1 :: (forall d. Data d => c (t d)) -> Maybe (c ChiSquared)
$cdataCast1 :: forall (t :: * -> *) (c :: * -> *).
Typeable t =>
(forall d. Data d => c (t d)) -> Maybe (c ChiSquared)
dataTypeOf :: ChiSquared -> DataType
$cdataTypeOf :: ChiSquared -> DataType
toConstr :: ChiSquared -> Constr
$ctoConstr :: ChiSquared -> Constr
gunfold :: (forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c ChiSquared
$cgunfold :: forall (c :: * -> *).
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c ChiSquared
gfoldl :: (forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> ChiSquared -> c ChiSquared
$cgfoldl :: forall (c :: * -> *).
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g) -> ChiSquared -> c ChiSquared
$cp1Data :: Typeable ChiSquared
Data, (forall x. ChiSquared -> Rep ChiSquared x)
-> (forall x. Rep ChiSquared x -> ChiSquared) -> Generic ChiSquared
forall x. Rep ChiSquared x -> ChiSquared
forall x. ChiSquared -> Rep ChiSquared x
forall a.
(forall x. a -> Rep a x) -> (forall x. Rep a x -> a) -> Generic a
$cto :: forall x. Rep ChiSquared x -> ChiSquared
$cfrom :: forall x. ChiSquared -> Rep ChiSquared x
Generic)

instance Show ChiSquared where
  showsPrec :: Int -> ChiSquared -> ShowS
showsPrec Int
i (ChiSquared Int
n) = String -> Int -> Int -> ShowS
forall a. Show a => String -> a -> Int -> ShowS
defaultShow1 String
"chiSquared" Int
n Int
i
instance Read ChiSquared where
  readPrec :: ReadPrec ChiSquared
readPrec = String -> (Int -> Maybe ChiSquared) -> ReadPrec ChiSquared
forall a r. Read a => String -> (a -> Maybe r) -> ReadPrec r
defaultReadPrecM1 String
"chiSquared" Int -> Maybe ChiSquared
chiSquaredE

instance ToJSON ChiSquared
instance FromJSON ChiSquared where
  parseJSON :: Value -> Parser ChiSquared
parseJSON (Object Object
v) = do
    Int
n <- Object
v Object -> Key -> Parser Int
forall a. FromJSON a => Object -> Key -> Parser a
.: Key
"chiSquaredNDF"
    Parser ChiSquared
-> (ChiSquared -> Parser ChiSquared)
-> Maybe ChiSquared
-> Parser ChiSquared
forall b a. b -> (a -> b) -> Maybe a -> b
maybe (String -> Parser ChiSquared
forall (m :: * -> *) a. MonadFail m => String -> m a
fail (String -> Parser ChiSquared) -> String -> Parser ChiSquared
forall a b. (a -> b) -> a -> b
$ Int -> String
errMsg Int
n) ChiSquared -> Parser ChiSquared
forall (m :: * -> *) a. Monad m => a -> m a
return (Maybe ChiSquared -> Parser ChiSquared)
-> Maybe ChiSquared -> Parser ChiSquared
forall a b. (a -> b) -> a -> b
$ Int -> Maybe ChiSquared
chiSquaredE Int
n
  parseJSON Value
_ = Parser ChiSquared
forall (f :: * -> *) a. Alternative f => f a
empty

instance Binary ChiSquared where
  put :: ChiSquared -> Put
put (ChiSquared Int
x) = Int -> Put
forall t. Binary t => t -> Put
put Int
x
  get :: Get ChiSquared
get = do Int
n <- Get Int
forall t. Binary t => Get t
get
           Get ChiSquared
-> (ChiSquared -> Get ChiSquared)
-> Maybe ChiSquared
-> Get ChiSquared
forall b a. b -> (a -> b) -> Maybe a -> b
maybe (String -> Get ChiSquared
forall (m :: * -> *) a. MonadFail m => String -> m a
fail (String -> Get ChiSquared) -> String -> Get ChiSquared
forall a b. (a -> b) -> a -> b
$ Int -> String
errMsg Int
n) ChiSquared -> Get ChiSquared
forall (m :: * -> *) a. Monad m => a -> m a
return (Maybe ChiSquared -> Get ChiSquared)
-> Maybe ChiSquared -> Get ChiSquared
forall a b. (a -> b) -> a -> b
$ Int -> Maybe ChiSquared
chiSquaredE Int
n


-- | Construct chi-squared distribution. Number of degrees of freedom
--   must be positive.
chiSquared :: Int -> ChiSquared
chiSquared :: Int -> ChiSquared
chiSquared Int
n = ChiSquared
-> (ChiSquared -> ChiSquared) -> Maybe ChiSquared -> ChiSquared
forall b a. b -> (a -> b) -> Maybe a -> b
maybe (String -> ChiSquared
forall a. HasCallStack => String -> a
error (String -> ChiSquared) -> String -> ChiSquared
forall a b. (a -> b) -> a -> b
$ Int -> String
errMsg Int
n) ChiSquared -> ChiSquared
forall a. a -> a
id (Maybe ChiSquared -> ChiSquared) -> Maybe ChiSquared -> ChiSquared
forall a b. (a -> b) -> a -> b
$ Int -> Maybe ChiSquared
chiSquaredE Int
n

-- | Construct chi-squared distribution. Number of degrees of freedom
--   must be positive.
chiSquaredE :: Int -> Maybe ChiSquared
chiSquaredE :: Int -> Maybe ChiSquared
chiSquaredE Int
n
  | Int
n Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
0    = Maybe ChiSquared
forall a. Maybe a
Nothing
  | Bool
otherwise = ChiSquared -> Maybe ChiSquared
forall a. a -> Maybe a
Just (Int -> ChiSquared
ChiSquared Int
n)

errMsg :: Int -> String
errMsg :: Int -> String
errMsg Int
n = String
"Statistics.Distribution.ChiSquared.chiSquared: N.D.F. must be positive. Got " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Int -> String
forall a. Show a => a -> String
show Int
n

instance D.Distribution ChiSquared where
  cumulative :: ChiSquared -> Double -> Double
cumulative = ChiSquared -> Double -> Double
cumulative

instance D.ContDistr ChiSquared where
  density :: ChiSquared -> Double -> Double
density ChiSquared
chi Double
x
    | Double
x Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
<= Double
0    = Double
0
    | Bool
otherwise = Double -> Double
forall a. Floating a => a -> a
exp (Double -> Double) -> Double -> Double
forall a b. (a -> b) -> a -> b
$ Double -> Double
forall a. Floating a => a -> a
log Double
x Double -> Double -> Double
forall a. Num a => a -> a -> a
* (Double
ndf2 Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
1) Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
x2 Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double -> Double
logGamma Double
ndf2 Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double -> Double
forall a. Floating a => a -> a
log Double
2 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
ndf2
    where
      ndf :: Double
ndf  = Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int -> Double) -> Int -> Double
forall a b. (a -> b) -> a -> b
$ ChiSquared -> Int
chiSquaredNDF ChiSquared
chi
      ndf2 :: Double
ndf2 = Double
ndfDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2
      x2 :: Double
x2   = Double
xDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2

  logDensity :: ChiSquared -> Double -> Double
logDensity ChiSquared
chi Double
x
    | Double
x Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
<= Double
0    = Double
m_neg_inf
    | Bool
otherwise = Double -> Double
forall a. Floating a => a -> a
log Double
x Double -> Double -> Double
forall a. Num a => a -> a -> a
* (Double
ndf2 Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
1) Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
x2 Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double -> Double
logGamma Double
ndf2 Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double -> Double
forall a. Floating a => a -> a
log Double
2 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
ndf2
    where
      ndf :: Double
ndf  = Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int -> Double) -> Int -> Double
forall a b. (a -> b) -> a -> b
$ ChiSquared -> Int
chiSquaredNDF ChiSquared
chi
      ndf2 :: Double
ndf2 = Double
ndfDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2
      x2 :: Double
x2   = Double
xDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2

  quantile :: ChiSquared -> Double -> Double
quantile = ChiSquared -> Double -> Double
quantile

instance D.Mean ChiSquared where
    mean :: ChiSquared -> Double
mean (ChiSquared Int
ndf) = Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
ndf

instance D.Variance ChiSquared where
    variance :: ChiSquared -> Double
variance (ChiSquared Int
ndf) = Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
ndf)

instance D.MaybeMean ChiSquared where
    maybeMean :: ChiSquared -> Maybe Double
maybeMean = Double -> Maybe Double
forall a. a -> Maybe a
Just (Double -> Maybe Double)
-> (ChiSquared -> Double) -> ChiSquared -> Maybe Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ChiSquared -> Double
forall d. Mean d => d -> Double
D.mean

instance D.MaybeVariance ChiSquared where
    maybeStdDev :: ChiSquared -> Maybe Double
maybeStdDev   = Double -> Maybe Double
forall a. a -> Maybe a
Just (Double -> Maybe Double)
-> (ChiSquared -> Double) -> ChiSquared -> Maybe Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ChiSquared -> Double
forall d. Variance d => d -> Double
D.stdDev
    maybeVariance :: ChiSquared -> Maybe Double
maybeVariance = Double -> Maybe Double
forall a. a -> Maybe a
Just (Double -> Maybe Double)
-> (ChiSquared -> Double) -> ChiSquared -> Maybe Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ChiSquared -> Double
forall d. Variance d => d -> Double
D.variance

instance D.Entropy ChiSquared where
  entropy :: ChiSquared -> Double
entropy (ChiSquared Int
ndf) =
    let kHalf :: Double
kHalf = Double
0.5 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
ndf in
    Double
kHalf
    Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double -> Double
forall a. Floating a => a -> a
log Double
2
    Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double -> Double
logGamma Double
kHalf
    Double -> Double -> Double
forall a. Num a => a -> a -> a
+ (Double
1Double -> Double -> Double
forall a. Num a => a -> a -> a
-Double
kHalf) Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double -> Double
digamma Double
kHalf

instance D.MaybeEntropy ChiSquared where
  maybeEntropy :: ChiSquared -> Maybe Double
maybeEntropy = Double -> Maybe Double
forall a. a -> Maybe a
Just (Double -> Maybe Double)
-> (ChiSquared -> Double) -> ChiSquared -> Maybe Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ChiSquared -> Double
forall d. Entropy d => d -> Double
D.entropy

instance D.ContGen ChiSquared where
    genContVar :: ChiSquared -> g -> m Double
genContVar (ChiSquared Int
n) = Int -> g -> m Double
forall g (m :: * -> *). StatefulGen g m => Int -> g -> m Double
MWC.chiSquare Int
n


cumulative :: ChiSquared -> Double -> Double
cumulative :: ChiSquared -> Double -> Double
cumulative ChiSquared
chi Double
x
  | Double
x Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
<= Double
0    = Double
0
  | Bool
otherwise = Double -> Double -> Double
incompleteGamma (Double
ndfDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2) (Double
xDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2)
  where
    ndf :: Double
ndf = Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int -> Double) -> Int -> Double
forall a b. (a -> b) -> a -> b
$ ChiSquared -> Int
chiSquaredNDF ChiSquared
chi

quantile :: ChiSquared -> Double -> Double
quantile :: ChiSquared -> Double -> Double
quantile (ChiSquared Int
ndf) Double
p
  | Double
p Double -> Double -> Bool
forall a. Eq a => a -> a -> Bool
== Double
0         = Double
0
  | Double
p Double -> Double -> Bool
forall a. Eq a => a -> a -> Bool
== Double
1         = Double
1Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
0
  | Double
p Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
> Double
0 Bool -> Bool -> Bool
&& Double
p Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
< Double
1 = Double
2 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double -> Double -> Double
invIncompleteGamma (Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
ndf Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Double
2) Double
p
  | Bool
otherwise      =
    String -> Double
forall a. HasCallStack => String -> a
error (String -> Double) -> String -> Double
forall a b. (a -> b) -> a -> b
$ String
"Statistics.Distribution.ChiSquared.quantile: p must be in [0,1] range. Got: "String -> ShowS
forall a. [a] -> [a] -> [a]
++Double -> String
forall a. Show a => a -> String
show Double
p