fallbacks.jl 7.01 KB
Newer Older
1 2 3 4 5 6
##############################################################################
#
# Fallback methods, usually overridden for specific distributions
#
##############################################################################

7 8 9 10 11 12 13
# array-like characteristics

# size(d::Distribution) = size(rand(d))
size(d::UnivariateDistribution) = ()
size(d::MultivariateDistribution) = (dim(d),)
size(d::MatrixDistribution) = (dim(d), dim(d)) # override if the matrix isn't square

Dahua Lin's avatar
Dahua Lin committed
14 15
length(d::UnivariateDistribution) = 1
length(d::MultivariateDistribution) = dim(d)
16 17
length(d::Distribution) = prod(size(d))

18 19 20 21 22 23
# support handling

isbounded(d::UnivariateDistribution) = isupperbounded(d) && islowerbounded(d)
hasfinitesupport(d::DiscreteUnivariateDistribution) = isbounded(d)
hasfinitesupport(d::ContinuousUnivariateDistribution) = false

24 25
insupport(d::Distribution, x) = false

26 27 28 29 30 31 32 33
function insupport!{D<:UnivariateDistribution}(r::AbstractArray, d::Union(D,Type{D}), X::AbstractArray)
    if length(r) != length(X)
        throw(ArgumentError("Inconsistent array dimensions."))
    end
    for i in 1 : length(X)
        r[i] = insupport(d, X[i])
    end
    r
34 35
end

36
insupport{D<:UnivariateDistribution}(d::Union(D,Type{D}), X::AbstractArray) = insupport!(BitArray(size(X)), d, X)
37

38 39 40 41 42 43 44
function insupport!{D<:MultivariateDistribution}(r::AbstractArray, d::Union(D,Type{D}), X::AbstractArray)
    n = div(length(X),size(X,1))
    if length(r) != n
        throw(ArgumentError("Inconsistent array dimensions."))
    end    
    for i in 1 : n
        r[i] = insupport(d, X[:, i])
45
    end
46
    r
47
end
48
insupport{D<:MultivariateDistribution}(d::Union(D,Type{D}), X::AbstractArray) = insupport!(BitArray(size(X)[2:end]), d, X)
49

50 51 52 53 54 55 56
function insupport!{D<:MatrixDistribution}(r::AbstractArray, d::Union(D,Type{D}), X::AbstractArray)
    n = div(length(X),size(X,1)*size(X,2))
    if length(r) != n
        throw(ArgumentError("Inconsistent array dimensions."))
    end    
    for i in 1 : n
        r[i] = insupport(d, X[:, :, i])
57
    end
58
    r
59
end
60
insupport{D<:MatrixDistribution}(d::Union(D,Type{D}), X::AbstractArray) = insupport!(BitArray(size(X)[3:end]), d, X)
61

62 63 64 65 66
# generic function to get number of samples

nsamples{D<:UnivariateDistribution}(dt::Type{D}, x::Array) = length(x)
nsamples{D<:MultivariateDistribution}(dt::Type{D}, x::Matrix) = size(x, 2)
nsamples{D<:MatrixDistribution,T}(dt::Type{D}, x::Array{Matrix{T}}) = length(x)
Dahua Lin's avatar
Dahua Lin committed
67 68

#### Statistics ####
69
mean(d::Distribution) = throw(MethodError(mean,(d,)))
Dahua Lin's avatar
Dahua Lin committed
70 71 72 73 74 75 76 77 78 79
std(d::Distribution) = sqrt(var(d))

# What's the purpose for this function?
function var{M <: Real}(d::UnivariateDistribution, mu::AbstractArray{M})
    V = similar(mu, Float64)
    for i in 1:length(mu)
        V[i] = var(d, mu[i])
    end
    return V
end
80

81
function cor(d::MultivariateDistribution)
82
    R = cov(d)
83 84 85 86 87 88 89 90 91
    m, n = size(R)
    for j in 1:n
        for i in 1:n
            R[i, j] = d.cov[i, j] / sqrt(d.cov[i, i] * d.cov[j, j])
        end
    end
    return R
end

Dahua Lin's avatar
Dahua Lin committed
92
binaryentropy(d::Distribution) = entropy(d) / log(2)
93

Dahua Lin's avatar
Dahua Lin committed
94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111
# kurtosis returns excess kurtosis by default
# proper kurtosis can be returned with correction = false
function kurtosis(d::Distribution, correction::Bool)
    if correction
        return kurtosis(d)
    else
        return kurtosis(d) + 3.0
    end
end

excess(d::Distribution) = kurtosis(d)
excess_kurtosis(d::Distribution) = kurtosis(d)
proper_kurtosis(d::Distribution) = kurtosis(d, false)

isplatykurtic(d::Distribution) = kurtosis(d) > 0.0
isleptokurtic(d::Distribution) = kurtosis(d) < 0.0
ismesokurtic(d::Distribution) = kurtosis(d) == 0.0

Dahua Lin's avatar
Dahua Lin committed
112
median(d::UnivariateDistribution) = quantile(d, 0.5)
113

114 115
modes(d::Distribution) = [mode(d)]

Dahua Lin's avatar
Dahua Lin committed
116 117 118 119 120 121
#### pdf, cdf, and quantile ####

logpdf(d::UnivariateDistribution, x::Real) = log(pdf(d, x))
pdf(d::MultivariateDistribution, x::Vector) = exp(logpdf(d, x))

ccdf(d::UnivariateDistribution, q::Real) = 1.0 - cdf(d, q)
122 123
cquantile(d::UnivariateDistribution, p::Real) = quantile(d, 1.0 - p)

Dahua Lin's avatar
Dahua Lin committed
124 125 126 127 128 129 130
logcdf(d::Distribution, q::Real) = log(cdf(d,q))
logccdf(d::Distribution, q::Real) = log(ccdf(d,q))
invlogccdf(d::Distribution, lp::Real) = quantile(d, -expm1(lp))
invlogcdf(d::Distribution, lp::Real) = quantile(d, exp(lp))


#### log likelihood ####
131

132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147
function loglikelihood(d::UnivariateDistribution, X::Array)
    ll = 0.0
    for i in 1:length(X)
        ll += logpdf(d, X[i])
    end
    return ll
end

function loglikelihood(d::MultivariateDistribution, X::Matrix)
    ll = 0.0
    for i in 1:size(X, 2)
        ll += logpdf(d, X[:, i])
    end
    return ll
end

148

Dahua Lin's avatar
Dahua Lin committed
149
#### Vectorized functions for univariate distributions ####
150

Dahua Lin's avatar
Dahua Lin committed
151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167
for fun in [:pdf, :logpdf, :cdf, :logcdf, :ccdf, :logccdf, :invlogcdf, :invlogccdf, :quantile, :cquantile]
    fun! = symbol(string(fun, '!'))

    @eval begin
        function ($fun!)(r::AbstractArray, d::UnivariateDistribution, X::AbstractArray)
            if length(r) != length(X)
                throw(ArgumentError("Inconsistent array dimensions."))
            end
            for i in 1 : length(X)
                r[i] = ($fun)(d, X[i])
            end
            r
        end

        function ($fun)(d::UnivariateDistribution, X::AbstractArray)
            ($fun!)(Array(Float64, size(X)), d, X)
        end
168 169 170
    end
end

Dahua Lin's avatar
Dahua Lin committed
171
#### Vectorized functions for multivariate distributions ####
172

173 174 175 176 177
function logpdf!(r::AbstractArray, d::MultivariateDistribution, x::AbstractMatrix)
    n::Int = size(x, 2)
    if length(r) != n
        throw(ArgumentError("Inconsistent array dimensions."))
    end
Dahua Lin's avatar
Dahua Lin committed
178
    for i in 1 : n
179 180
        r[i] = logpdf(d, x[:, i])
    end
Dahua Lin's avatar
Dahua Lin committed
181
    r
182 183
end

Dahua Lin's avatar
Dahua Lin committed
184 185
function logpdf(d::MultivariateDistribution, X::AbstractMatrix)  
    logpdf!(Array(Float64, size(X, 2)), d, X)
186 187
end

Dahua Lin's avatar
Dahua Lin committed
188 189 190 191 192 193 194
function pdf!(r::AbstractArray, d::MultivariateDistribution, X::AbstractMatrix)
    logpdf!(r, d, X)  # size checking is done by logpdf!
    for i in 1 : size(X, 2)
        r[i] = exp(r[i])
    end
    r
end
195

196
function pdf(d::MultivariateDistribution, X::AbstractMatrix)
Dahua Lin's avatar
Dahua Lin committed
197
    pdf!(Array(Float64, size(X, 2)), d, X)
198 199 200
end


Dahua Lin's avatar
Dahua Lin committed
201 202 203 204 205 206 207 208
#### logpmf & pmf for discrete distributions ####

logpmf(d::DiscreteDistribution, args::Any...) = logpdf(d, args...)
logpmf!(r::AbstractArray, d::DiscreteDistribution, args::Any...) = logpdf!(r, d, args...)
pmf(d::DiscreteDistribution, args::Any...) = pdf(d, args...)


#### Sampling: rand & rand! ####
209

210 211 212
# default: inverse transform sampling
rand(d::UnivariateDistribution) = quantile(d, rand())

213 214 215 216
function sprand(m::Integer, n::Integer, density::Real, d::Distribution)
    return sprand(m, n, density, n -> rand(d, n))
end

Dahua Lin's avatar
Dahua Lin committed
217

218 219
# Fitting

220 221 222 223 224
function suffstats{D<:Distribution}(dt::Type{D}, xs...) 
    argtypes = tuple(D, map(typeof, xs)...)
    error("suffstats is not implemented for $argtypes.")
end

Dahua Lin's avatar
Dahua Lin committed
225 226 227 228 229 230
fit_mle{D<:UnivariateDistribution}(dt::Type{D}, x::Array) = fit_mle(D, suffstats(D, x))
fit_mle{D<:UnivariateDistribution}(dt::Type{D}, x::Array, w::Array) = fit_mle(D, suffstats(D, x, w))

fit_mle{D<:MultivariateDistribution}(dt::Type{D}, x::Matrix) = fit_mle(D, suffstats(D, x))
fit_mle{D<:MultivariateDistribution}(dt::Type{D}, x::Matrix, w::Array) = fit_mle(D, suffstats(D, x, w))

231 232
fit{D<:Distribution}(dt::Type{D}, x) = fit_mle(D, x)
fit{D<:Distribution}(dt::Type{D}, args...) = fit_mle(D, args...)
233