Commit 244c27c0 authored by Simon Byrne's avatar Simon Byrne

fix sqrt constants now √ an operator

parent d835cb66
......@@ -3,8 +3,8 @@
import Base.@math_const
@math_const twoπ 6.2831853071795864769 big(2.) * π
@math_const 2 1.4142135623730950488 sqrt(big(2.))
@math_const sqrt2 1.4142135623730950488 sqrt(big(2.))
@math_const log2π 1.8378770664093454836 log(big(2.)*π)
@math_const 2.5066282746310005024 sqrt(big(2.)*π)
@math_const sqrt 2.5066282746310005024 sqrt(big(2.)*π)
@math_const logtwo 0.6931471805599453094 log(big(2.))
@math_const loghalf -0.6931471805599453094 log(big(0.5))
......@@ -21,13 +21,13 @@ logexpm1(x::Integer) = logexpm1(float(x))
# log(1+x^2)
log1psq(x::FloatingPoint) = (ax = abs(x); ax < maxintfloat(x) ? log1p(ax*ax) : 2*log(ax))
φ(z::Real) = exp(-0.5*z*z)/
φ(z::Real) = exp(-0.5*z*z)/sqrt
logφ(z::Real) = -0.5*(z*z + log2π)
Φ(z::Real) = 0.5*erfc(-z/2)
Φc(z::Real) = 0.5*erfc(z/2)
logΦ(z::Real) = z < -1.0 ? log(0.5*erfcx(-z/√2)) - 0.5*z*z : log1p(-0.5*erfc(z/2))
logΦc(z::Real) = z > 1.0 ? log(0.5*erfcx(z/√2)) - 0.5*z*z : log1p(-0.5*erfc(-z/2))
Φ(z::Real) = 0.5*erfc(-z/sqrt2)
Φc(z::Real) = 0.5*erfc(z/sqrt2)
logΦ(z::Real) = z < -1.0 ? log(0.5*erfcx(-z/sqrt2)) - 0.5*z*z : log1p(-0.5*erfc(z/sqrt2))
logΦc(z::Real) = z > 1.0 ? log(0.5*erfcx(z/sqrt2)) - 0.5*z*z : log1p(-0.5*erfc(-z/sqrt2))
import Base.Math.@horner
......
......@@ -18,7 +18,7 @@ cquantile(d::Chi,p::Real) = sqrt(cquantile(Chisq(d.df),p))
invlogcdf(d::Chi,p::Real) = sqrt(invlogcdf(Chisq(d.df),p))
invlogccdf(d::Chi,p::Real) = sqrt(invlogccdf(Chisq(d.df),p))
mean(d::Chi) = 2 * gamma((d.df + 1.0) / 2.0) / gamma(d.df / 2.0)
mean(d::Chi) = sqrt2 * gamma((d.df + 1.0) / 2.0) / gamma(d.df / 2.0)
function mode(d::Chi)
d.df >= 1.0 || error("Chi distribution has no mode when df < 1")
......
......@@ -7,7 +7,7 @@ end
@continuous_distr_support Kolmogorov 0.0 Inf
mean(d::Kolmogorov) = 0.5**log(2.0)
mean(d::Kolmogorov) = 0.5*sqrt*log(2.0)
var(d::Kolmogorov) = pi*pi/12.0 - 0.5*pi*log(2.0)^2
# TODO: higher-order moments also exist, can be obtained by differentiating series
......@@ -31,7 +31,7 @@ function cdf(d::Kolmogorov,x::Real)
for i = 1:5
s += exp(-((2*i-1)*π/x)^2/8.0)
end
*s/x
sqrt*s/x
end
function ccdf(d::Kolmogorov,x::Real)
......@@ -56,7 +56,7 @@ function pdf(d::Kolmogorov,x::Real)
k = ((2*i-1)*c)^2
s += (k-1.0)*exp(-k/2.0)
end
return *s/x^2
return sqrt*s/x^2
else
s = 0.0
for i = 1:20
......
......@@ -46,7 +46,7 @@ end
function pdf(d::Levy, x::Real)
m, c = d.location, d.scale
(sqrt(c)/) * exp(-(c / (2.0 * (x - m)))) / (x - m)^1.5
(sqrt(c)/sqrt) * exp(-(c / (2.0 * (x - m)))) / (x - m)^1.5
end
quantile(d::Levy, p::Real) = d.location + d.scale / (2.erfcinv(p)^2)
......
......@@ -9,7 +9,7 @@ end
mean(d::NoncentralChisq) = d.df + d.ncp
var(d::NoncentralChisq) = 2.0*(d.df + 2.0*d.ncp)
skewness(d::NoncentralChisq) = 2.0*2*(d.df + 3.0*d.ncp)/sqrt(d.df + 2.0*d.ncp)^3
skewness(d::NoncentralChisq) = 2.0*sqrt2*(d.df + 3.0*d.ncp)/sqrt(d.df + 2.0*d.ncp)^3
kurtosis(d::NoncentralChisq) = 12.0*(d.df + 4.0*d.ncp)/(d.df + 2.0*d.ncp)^2
function mgf(d::NoncentralChisq, t::Real)
......
......@@ -38,7 +38,7 @@ entropy(cf::NormalCanon) = 0.5 * (log2π + 1. - log(cf.prec))
# Evaluation
pdf(d::NormalCanon, x::Real) = (sqrt(d.prec) / ) * exp(-0.5 * d.prec * abs2(x - d.μ))
pdf(d::NormalCanon, x::Real) = (sqrt(d.prec) / sqrt) * exp(-0.5 * d.prec * abs2(x - d.μ))
logpdf(d::NormalCanon, x::Real) = 0.5 * (log(d.prec) - log2π - d.prec * abs2(x - d.μ))
zval(d::NormalCanon, x::Real) = (x - d.μ) * sqrt(d.prec)
......
......@@ -52,7 +52,7 @@ end
function randnt(lower::Real, upper::Real)
if (lower <= 0 && upper == Inf) ||
(upper >= 0 && lower == Inf) ||
(lower <= 0 && upper >= 0 && upper - lower > )
(lower <= 0 && upper >= 0 && upper - lower > sqrt)
while true
r = randn()
if r > lower && r < upper
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment