julia> using WildBootTests, CSV, DataFrames, StatsModels, GLM, Plots
julia> d = download("https://raw.github.com/vincentarelbundock/Rdatasets/master/csv/sandwich/PetersenCL.csv");
julia> df = CSV.read(d, DataFrame);
julia> f = @formula(y ~ 1 + x); # state OLS model
julia> f = apply_schema(f, schema(f, df)); # link model to data
julia> lm(f, df) # run OLS for illustration; not needed for following lines
StatsModels.TableRegressionModel{LinearModel{GLM.LmResp{Vector{Float64}}, GLM.DensePredChol{Float64, LinearAlgebra.CholeskyPivoted{Float64, Matrix{Float64}}}}, Matrix{Float64}}
y ~ 1 + x
Coefficients:
─────────────────────────────────────────────────────────────────────────
Coef. Std. Error t Pr(>|t|) Lower 95% Upper 95%
─────────────────────────────────────────────────────────────────────────
(Intercept) 0.0296797 0.0283593 1.05 0.2954 -0.025917 0.0852764
x 1.03483 0.0285833 36.20 <1e-99 0.978798 1.09087
─────────────────────────────────────────────────────────────────────────
julia> resp, predexog = modelcols(f, df); # extract response & (exogenous) predictor variables
julia> clustid = df.firm; # extract clustering variable
julia> R = [0 1]; r = [1]; # put null that coefficient on x = 1 in Rβ̂ = r form, where β̂ is parameter vector
julia> test = wildboottest(R, r; resp=resp, predexog=predexog, clustid=clustid)
WildBootTests.BootTestResult{Float32}
p = 0.492
ci = Float32[0.93461335 1.1347668]
julia> test = wildboottest(R, r; resp, predexog, clustid); # same, using Julia syntactic sugar
julia> p(test) # programmatically extract p value
0.49459493f0
julia> ci(test) # programmatically extract confidence interval
1×2 Matrix{Float32}:
0.934961 1.13469
julia> plot(plotpoints(test)...) # plot confidence curve
using WildBootTests, CSV, DataFrames, StatsModels, GLM, Plots
# use Webb instead of Rademacher weights, 99,999 bootstrap replications instead of 999
wildboottest(R, r; resp, predexog, clustid, reps=99999, auxwttype=:webb)
# jackknife the bootstrap data-generating process to reduce distortion from outliers
wildboottest(Float64, R, r; resp, predexog, clustid, jk=true)
# use guaranteed-stable random number generator for exact replicability
using StableRNGs
wildboottest(R, r; resp, predexog, clustid, rng=StableRNG(23948572))
# test that coefficient on intercept = 0 and coefficient on x = 1; plot confidence surface
test = wildboottest([1 0; 0 1], [0;1]; resp, predexog, clustid, reps=9999)
plot(plotpoints(test).X..., plotpoints(test).p, st=:contourf)
# multiway-cluster errors by firm and year; bootstrap by firm
wildboottest(R, r; resp, predexog, clustid=Matrix(df[:,[:firm, :year]]), nerrclustvar=2, nbootclustvar=1)
# same but bootstrap by year
wildboottest(R, r; resp, predexog, clustid=Matrix(df[:,[:year, :firm]]), nerrclustvar=2, nbootclustvar=1)
# same but bootstrap by year-firm pair
wildboottest(R, r; resp, predexog, clustid=Matrix(df[:,[:year, :firm]]), nerrclustvar=2, nbootclustvar=2)
# Rao/score test with multiway clustering of errors but no bootstrap
wildboottest(R, r; resp, predexog, predendog, inst, Matrix(df[:,[:year, :firm]]), reps=0)
# Same but Wald test: i.e., conventional, multiway clustered errors
wildboottest(R, r; resp, predexog, predendog, inst, clustid=Matrix(df[:,[:year, :firm]]), reps=0, imposenull=false)
# add year fixed effects to model; cluster by firm
wildboottest(R, r; resp, predexog, feid=df.year, clustid=df.firm)
# test hypotheses, while imposing model constraint that constant term = 0.2
R1 = [1 0]; r1 = [.2]
wildboottest(R, r; R1, r1, resp, predexog, clustid=df.firm)