-- |-- Module : Statistics.Resampling.Bootstrap-- Copyright : (c) 2009, 2011 Bryan O'Sullivan-- License : BSD3---- Maintainer : bos@serpentine.com-- Stability : experimental-- Portability : portable---- The bootstrap method for statistical inference.moduleStatistics.Resampling.Bootstrap(bootstrapBCA ,basicBootstrap -- * References-- $references)whereimportData.Vector.Generic((!))importqualifiedData.Vector.UnboxedasUimportqualifiedData.Vector.GenericasGimportStatistics.Distribution (cumulative ,quantile )importStatistics.Distribution.Normal importStatistics.Resampling (Bootstrap (..),jackknife )importStatistics.Sample (mean )importStatistics.Types (Sample ,CL ,Estimate ,ConfInt ,estimateFromInterval ,estimateFromErr ,CL ,significanceLevel )importStatistics.Function (gsort )importqualifiedStatistics.Resampling asRimportControl.Parallel.Strategies(parMap,rdeepseq)dataT ={-# UNPACK#-}!Double:< {-# UNPACK#-}!Doubleinfixl2:< -- | Bias-corrected accelerated (BCA) bootstrap. This adjusts for both-- bias and skewness in the resampled distribution.---- BCA algorithm is described in ch. 5 of Davison, Hinkley "Confidence-- intervals" in section 5.3 "Percentile method"bootstrapBCA ::CL Double-- ^ Confidence level->Sample -- ^ Full data sample->[(R.Estimator ,Bootstrap U.VectorDouble)]-- ^ Estimates obtained from resampled data and estimator used for-- this.->[Estimate ConfInt Double]bootstrapBCA :: CL Double
-> Sample
-> [(Estimator, Bootstrap Vector Double)]
-> [Estimate ConfInt Double]
bootstrapBCA CL Double
confidenceLevel Sample
sample [(Estimator, Bootstrap Vector Double)]
resampledData =Strategy (Estimate ConfInt Double)
-> ((Estimator, Bootstrap Vector Double)
 -> Estimate ConfInt Double)
-> [(Estimator, Bootstrap Vector Double)]
-> [Estimate ConfInt Double]
forall b a. Strategy b -> (a -> b) -> [a] -> [b]
parMapStrategy (Estimate ConfInt Double)
forall a. NFData a => Strategy a
rdeepseq(Estimator, Bootstrap Vector Double) -> Estimate ConfInt Double
forall {a}.
(Num a, Unbox a, Ord a) =>
(Estimator, Bootstrap Vector a) -> Estimate ConfInt a
e [(Estimator, Bootstrap Vector Double)]
resampledData wheree :: (Estimator, Bootstrap Vector a) -> Estimate ConfInt a
e (Estimator
est ,Bootstrap a
pt Vector a
resample )|Sample -> Int
forall a. Unbox a => Vector a -> Int
U.lengthSample
sample Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
==Int
1Bool -> Bool -> Bool
||Double -> Bool
forall a. RealFloat a => a -> Bool
isInfiniteDouble
bias =a -> (a, a) -> CL Double -> Estimate ConfInt a
forall a. a -> (a, a) -> CL Double -> Estimate ConfInt a
estimateFromErr a
pt (a
0,a
0)CL Double
confidenceLevel |Bool
otherwise=a -> (a, a) -> CL Double -> Estimate ConfInt a
forall a. Num a => a -> (a, a) -> CL Double -> Estimate ConfInt a
estimateFromInterval a
pt (Vector a
resample Vector a -> Int -> a
forall (v :: * -> *) a.
(HasCallStack, Vector v a) =>
v a -> Int -> a
!Int
lo ,Vector a
resample Vector a -> Int -> a
forall (v :: * -> *) a.
(HasCallStack, Vector v a) =>
v a -> Int -> a
!Int
hi )CL Double
confidenceLevel where-- Quantile estimates for given CLlo :: Int
lo =Int -> Int -> Int
forall a. Ord a => a -> a -> a
min(Int -> Int -> Int
forall a. Ord a => a -> a -> a
max(Double -> Int
cumn Double
a1 )Int
0)(Int
ni Int -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1)wherea1 :: Double
a1 =Double
bias Double -> Double -> Double
forall a. Num a => a -> a -> a
+Double
b1 Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/(Double
1Double -> Double -> Double
forall a. Num a => a -> a -> a
-Double
accel Double -> Double -> Double
forall a. Num a => a -> a -> a
*Double
b1 )b1 :: Double
b1 =Double
bias Double -> Double -> Double
forall a. Num a => a -> a -> a
+Double
z1 hi :: Int
hi =Int -> Int -> Int
forall a. Ord a => a -> a -> a
max(Int -> Int -> Int
forall a. Ord a => a -> a -> a
min(Double -> Int
cumn Double
a2 )(Int
ni Int -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1))Int
0wherea2 :: Double
a2 =Double
bias Double -> Double -> Double
forall a. Num a => a -> a -> a
+Double
b2 Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/(Double
1Double -> Double -> Double
forall a. Num a => a -> a -> a
-Double
accel Double -> Double -> Double
forall a. Num a => a -> a -> a
*Double
b2 )b2 :: Double
b2 =Double
bias Double -> Double -> Double
forall a. Num a => a -> a -> a
-Double
z1 -- Number of resamplesni :: Int
ni =Vector a -> Int
forall a. Unbox a => Vector a -> Int
U.lengthVector a
resample n :: Double
n =Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegralInt
ni -- Correctionsz1 :: Double
z1 =NormalDistribution -> Double -> Double
forall d. ContDistr d => d -> Double -> Double
quantile NormalDistribution
standard (CL Double -> Double
forall a. CL a -> a
significanceLevel CL Double
confidenceLevel Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2)cumn :: Double -> Int
cumn =Double -> Int
forall b. Integral b => Double -> b
forall a b. (RealFrac a, Integral b) => a -> b
round(Double -> Int) -> (Double -> Double) -> Double -> Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
.(Double -> Double -> Double
forall a. Num a => a -> a -> a
*Double
n )(Double -> Double) -> (Double -> Double) -> Double -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
.NormalDistribution -> Double -> Double
forall d. Distribution d => d -> Double -> Double
cumulative NormalDistribution
standard bias :: Double
bias =NormalDistribution -> Double -> Double
forall d. ContDistr d => d -> Double -> Double
quantile NormalDistribution
standard (Double
probN Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
n )whereprobN :: Double
probN =Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral(Int -> Double) -> (Vector a -> Int) -> Vector a -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
.Vector a -> Int
forall a. Unbox a => Vector a -> Int
U.length(Vector a -> Int) -> (Vector a -> Vector a) -> Vector a -> Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
.(a -> Bool) -> Vector a -> Vector a
forall a. Unbox a => (a -> Bool) -> Vector a -> Vector a
U.filter(a -> a -> Bool
forall a. Ord a => a -> a -> Bool
<a
pt )(Vector a -> Double) -> Vector a -> Double
forall a b. (a -> b) -> a -> b
$Vector a
resample accel :: Double
accel =Double
sumCubes Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/(Double
6Double -> Double -> Double
forall a. Num a => a -> a -> a
*(Double
sumSquares Double -> Double -> Double
forall a. Floating a => a -> a -> a
**Double
1.5))where(Double
sumSquares :< Double
sumCubes )=(T -> Double -> T) -> T -> Sample -> T
forall b a. Unbox b => (a -> b -> a) -> a -> Vector b -> a
U.foldl'T -> Double -> T
f (Double
0Double -> Double -> T
:< Double
0)Sample
jack f :: T -> Double -> T
f (Double
s :< Double
c )Double
j =Double
s Double -> Double -> Double
forall a. Num a => a -> a -> a
+Double
d2 Double -> Double -> T
:< Double
c Double -> Double -> Double
forall a. Num a => a -> a -> a
+Double
d2 Double -> Double -> Double
forall a. Num a => a -> a -> a
*Double
d whered :: Double
d =Double
jackMean Double -> Double -> Double
forall a. Num a => a -> a -> a
-Double
j d2 :: Double
d2 =Double
d Double -> Double -> Double
forall a. Num a => a -> a -> a
*Double
d jackMean :: Double
jackMean =Sample -> Double
forall (v :: * -> *). Vector v Double => v Double -> Double
mean Sample
jack jack :: Sample
jack =Estimator -> Sample -> Sample
jackknife Estimator
est Sample
sample -- | Basic bootstrap. This method simply uses empirical quantiles for-- confidence interval.basicBootstrap ::(G.Vectorv a ,Orda ,Numa )=>CL Double-- ^ Confidence vector->Bootstrap v a -- ^ Estimate from full sample and vector of-- estimates obtained from resamples->Estimate ConfInt a {-# INLINEbasicBootstrap #-}basicBootstrap :: forall (v :: * -> *) a.
(Vector v a, Ord a, Num a) =>
CL Double -> Bootstrap v a -> Estimate ConfInt a
basicBootstrap CL Double
cl (Bootstrap a
e v a
ests )=a -> (a, a) -> CL Double -> Estimate ConfInt a
forall a. Num a => a -> (a, a) -> CL Double -> Estimate ConfInt a
estimateFromInterval a
e (v a
sorted v a -> Int -> a
forall (v :: * -> *) a.
(HasCallStack, Vector v a) =>
v a -> Int -> a
!Int
lo ,v a
sorted v a -> Int -> a
forall (v :: * -> *) a.
(HasCallStack, Vector v a) =>
v a -> Int -> a
!Int
hi )CL Double
cl wheresorted :: v a
sorted =v a -> v a
forall e (v :: * -> *). (Ord e, Vector v e) => v e -> v e
gsort v a
ests n :: Double
n =Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral(Int -> Double) -> Int -> Double
forall a b. (a -> b) -> a -> b
$v a -> Int
forall (v :: * -> *) a. Vector v a => v a -> Int
G.lengthv a
ests c :: Double
c =Double
n Double -> Double -> Double
forall a. Num a => a -> a -> a
*(CL Double -> Double
forall a. CL a -> a
significanceLevel CL Double
cl Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2)-- FIXME: can we have better estimates of quantiles in case when p-- is not multiple of 1/N---- FIXME: we could have undercoverage herelo :: Int
lo =Double -> Int
forall b. Integral b => Double -> b
forall a b. (RealFrac a, Integral b) => a -> b
roundDouble
c hi :: Int
hi =Double -> Int
forall b. Integral b => Double -> b
forall a b. (RealFrac a, Integral b) => a -> b
truncate(Double
n Double -> Double -> Double
forall a. Num a => a -> a -> a
-Double
c )-- $references---- * Davison, A.C; Hinkley, D.V. (1997) Bootstrap methods and their-- application. <http://statwww.epfl.ch/davison/BMA/>