t_test.jl

using Distributions
using HypothesisTests

n_sims = 1_000_000
n = 20
x = Array{Float64}(n)
y = Array{Float64}(n)
p_d = Array{Float64}(n_sims)
for s in 1:n_sims
    for i in 1:n
        z = rand(Bernoulli(0.5))
        x[i] = z * (1e-12 * rand(Bernoulli(0.5))) + (1 - z) * rand(Bernoulli(1e-6))
        y[i] = z * (1e-12 * rand(Bernoulli(0.5))) + (1 - z) * rand(Bernoulli(1e-6))
    end
    try
        p_d[s] = pvalue(UnequalVarianceTTest(x, y))
    catch
        p_d[s] = NaN
    end
end
mean(isfinite.(p_d))
mean(p_d .< 0.05)
mean(p_d .< 0.01)
mean(p_d .< 0.001)

julia> mean(isfinite.(p_d))
0.999917

julia> mean(p_d .< 0.05)
0.018518

julia> mean(p_d .< 0.01)
0.001566

julia> mean(p_d .< 0.001)
4.9e-5

Bug here where z is re-used.