#sample size, number of covariates, number of true causative predictors, Distribution, Link
n = 1000
p = 10000
k = 10
d = Poisson
l = LogLink()

#set random seed
Random.seed!(2019)

#random normal design matrix
X = randn(n, p)

#define true model with non-zero entries 0.1, 0.2, ..., 1.0
true_b = zeros(p)
true_b[1:10] .= collect(0.1:0.1:1.0)
shuffle!(true_b)
correct_position = findall(!iszero, true_b)

#simulate y
μ = X * true_b
prob = linkinv.(l, μ)
clamp!(prob, -20, 20)
y = [rand(d(i)) for i in prob]
y = Float64.(y)

#now perform univariate tests to test for association
data = DataFrame(X=zeros(n), Y=y)
placeholder = zeros(n)
pvalues = zeros(p)
for i in 1:p
    copyto!(placeholder, @view(X[:, i]))
    data.X .= placeholder
    result = glm(@formula(Y ~ X), data, Poisson(), LogLink())
    pvalues[i] = coeftable(result).cols[4][2].v
end

#significance = bonferonni correction
significance = 0.05 / p
passing_position = findall(pvalues .<= significance)

#determine whether the 10 causative predictors were found
marginal_found = [correct_position[i] in passing_position for i in 1:k]

#compute some summary statistic
true_positive = count(!iszero, true_b[passing_snps_position2])
false_positive = length(passing_snps_position2) - juliaglm_tp
false_negative = k - juliaglm_tp

@show true_positive
@show false_positive
@show false_negative