Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 6 additions & 0 deletions .gitignore
Original file line number Diff line number Diff line change
Expand Up @@ -7,3 +7,9 @@
.vscode/
*~
*jld
*.jld2
/exper/Manifest*.toml
/exper/sim*/results/
/exper/sim*/*.csv
/exper/sim*/*.tex
/exper/sim*/*.pdf
4 changes: 2 additions & 2 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -44,8 +44,8 @@ CI = simulCIstudentized(straps, 0.95)
`Distributions` also exports `estimate`; when both packages are loaded, call
`DiffPrecTest.estimate`.

The `exper/` folder holds simulation scripts written for Julia 1.5 and version
0.1 of the package; they are kept as a record.
The `exper/` folder holds the simulation scripts, with their own environment;
see [`exper/README.md`](exper/README.md).

Package depends on [ProximalBase.jl](https://github.com/mlakolar/ProximalBase.jl),
[CoordinateDescent.jl](https://github.com/mlakolar/CoordinateDescent.jl) and
Expand Down
33 changes: 33 additions & 0 deletions exper/Project.toml
Original file line number Diff line number Diff line change
@@ -0,0 +1,33 @@
[deps]
CSV = "336ed68f-0bac-5ca0-87d4-7b16caf5d00b"
CovSel = "f50c8ead-2a5d-5d7a-afef-a6d87c1f6489"
DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
DiffPrecTest = "bb4bfe13-bfd7-5c6f-a2d2-8072376c6806"
Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
JLD2 = "033835bb-8acc-5ee8-8aae-3f567f8a3819"
LaTeXStrings = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f"
LaTeXTabulars = "266f59ce-6e72-579c-98bb-27b39b5c037e"
LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
ProximalBase = "0a755509-c5b1-5d6f-9b24-59bcda3c37a3"
Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"

[sources]
CovSel = {url = "https://github.com/mlakolar/CovSel.jl", rev = "master"}
DiffPrecTest = {path = ".."}

[compat]
CSV = "1"
CovSel = "0.2"
DataFrames = "1"
DiffPrecTest = "0.2"
Distributions = "0.25"
JLD2 = "0.6"
LaTeXStrings = "1"
LaTeXTabulars = "1"
Plots = "1"
ProximalBase = "0.3, 0.4"
StatsBase = "0.34"
julia = "1.10"
67 changes: 67 additions & 0 deletions exper/README.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,67 @@
# Simulations

Simulation scripts for DiffPrecTest. Each `simN` folder holds the scripts of one
simulation; one run of a script computes one repetition, and the Slurm jobs run
1000 repetitions (100 for sim5).

| Simulation | Question | Model | Scripts |
|---|---|---|---|
| sim1, sim2 | coverage of the confidence interval for one entry, `Δ = 0` | `Ω` banded (sim1), `Ω[i, j] = 0.9^\|i-j\|` (sim2) | `simN.jl`, `summarize.jl`, `plot_results.jl` (sim1) |
| sim3, sim4 | coverage for one entry, sparse `Δ` | `Ωx` banded (sim3), `Ωx[i, j] = 0.7^\|i-j\|` (sim4) | `simN.jl`, `summarize.jl` |
| sim6, sim7 | coverage for one entry, banded `Δ` | banded `Ωx` and `Ωy` | `simN.jl`, `summarize.jl` |
| sim5 | support recovery | block diagonal `Ωx`, diagonal `Ωy` (`model.jl`), `p = 60` (a) and `100` (b) | `sim5a.jl`, `sim5b.jl`, `summarize.jl` |
| sim8–sim11 | simultaneous bootstrap intervals for all entries, `Δ = 0` | banded, hub, random sparse and block `Ωx` | `simN.jl` (estimated supports), `simN_o.jl` (oracle supports), `simN_YinXia.jl` |
| sim8–sim12 | power of the bootstrap test along `Ωy = (1 - α) Ωx + α Ωy₁` | as above; banded `Ωx`, `Ωy` (sim12) | `simN_o_power.jl`, `simN_YinXia_power.jl` |

The coverage simulations compare nine methods: the three stage procedure
(`Sym-N`, `Asym-N`), Xia, Cai and Cai (2015) (`YinXia`), bootstrap standard
errors on the estimated support (`Sym-B`, `Asym-B`), and the same procedures on
the true support (`O-...`).

## Setup

With Julia 1.11 or later:

```sh
julia --project=. -e 'using Pkg; Pkg.instantiate()'
```

`Project.toml` takes DiffPrecTest from this repository and CovSel from GitHub.
Julia 1.10 ignores these sources; there, run
`using Pkg; Pkg.develop(path=".."); Pkg.add(url="https://github.com/mlakolar/CovSel.jl")`
once instead.

## Running

One repetition, from a simulation folder:

```sh
cd sim1
julia --project=.. sim1.jl 1 1 1 /tmp/sim1 # rep ip iElem dir
```

All repetitions on a Slurm cluster, with GNU parallel running one repetition per
task (`./start_all.sh` submits sims 1-7):

```sh
export SCRATCH=/path/to/scratch JULIA=/path/to/julia
cd sim1
./start_sbatch.sh -p <partition> # arguments go to sbatch
```

Results go to `$SCRATCH/diffTest/<simulation>`. The `simN_YinXia*.jl` scripts
run all repetitions in one process and write to the current folder; set
`NUM_REP` for fewer.

## Summaries

```sh
./collect.sh # sims 1-4, 6, 7 into sim*/results
julia --project=. sim1/summarize.jl # sim1.csv, sim1_table.tex
julia --project=. sim1/plot_results.jl # Q-Q plots and histograms
julia --project=. sim5/summarize.jl $SCRATCH/diffTest/sim5a a
julia --project=. sim8/collect_our.jl # bootstrap coverage, sims 8-11
```

The scripts were first run with Julia 1.5 in 2019–2020. Julia's default random
number generator changed in Julia 1.7, so the same seeds give different data.
24 changes: 5 additions & 19 deletions exper/cleanup.sh
100644 → 100755
Original file line number Diff line number Diff line change
@@ -1,21 +1,7 @@
#!/bin/bash

rm *~
rm sim1/*log
rm sim1/*out
rm sim1/*~
rm sim2/*log
rm sim2/*out
rm sim2/*~
rm sim3/*log
rm sim3/*out
rm sim3/*~
rm sim4/*log
rm sim4/*out
rm sim4/*~
# rm sim5/*log
# rm sim5/*out
# rm sim5/*~
rm sim6/*log
rm sim6/*out
rm sim6/*~
## removes the GNU parallel job logs and editor backups written next to the
## scripts

cd "$(dirname "$0")"
rm -f *~ sim*/*~ sim*/__*.runtask.log
120 changes: 33 additions & 87 deletions exper/collect.sh
100644 → 100755
Original file line number Diff line number Diff line change
@@ -1,87 +1,33 @@
#sim 6

# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim6/res_1_1 sim6_res_1_1
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim6/res_1_2 sim6_res_1_2
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim6/res_1_3 sim6_res_1_3
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim6/res_1_4 sim6_res_1_4
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim6/res_1_5 sim6_res_1_5

# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim6/res_2_1 sim6_res_2_1
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim6/res_2_2 sim6_res_2_2
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim6/res_2_3 sim6_res_2_3
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim6/res_2_4 sim6_res_2_4
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim6/res_2_5 sim6_res_2_5






# sim1
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim1/res_1_1 sim1_res_1_1
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim1/res_1_2 sim1_res_1_2
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim1/res_1_3 sim1_res_1_3
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim1/res_2_1 sim1_res_2_1
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim1/res_2_2 sim1_res_2_2
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim1/res_2_3 sim1_res_2_3


# sim2
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim2/res_1_1 sim2_res_1_1
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim2/res_1_2 sim2_res_1_2
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim2/res_1_3 sim2_res_1_3
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim2/res_2_1 sim2_res_2_1
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim2/res_2_2 sim2_res_2_2
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim2/res_2_3 sim2_res_2_3

# sim3
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_1_1 sim3_res_1_1
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_1_2 sim3_res_1_2
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_1_3 sim3_res_1_3
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_1_4 sim3_res_1_4
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_1_5 sim3_res_1_5

# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_2_1 sim3_res_2_1
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_2_2 sim3_res_2_2
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_2_3 sim3_res_2_3
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_2_4 sim3_res_2_4
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_2_5 sim3_res_2_5

# # jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_3_1 sim3_res_3_1
# # jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_3_2 sim3_res_3_2
# # jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_3_3 sim3_res_3_3
# # jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_3_4 sim3_res_3_4
# # jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim3/res_3_5 sim3_res_3_5


# sim4
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim4/res_1_1 sim4_res_1_1
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim4/res_1_2 sim4_res_1_2
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim4/res_1_3 sim4_res_1_3
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim4/res_1_4 sim4_res_1_4
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim4/res_1_5 sim4_res_1_5

# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim4/res_2_1 sim4_res_2_1
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim4/res_2_2 sim4_res_2_2
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim4/res_2_3 sim4_res_2_3
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim4/res_2_4 sim4_res_2_4
# jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim4/res_2_5 sim4_res_2_5





# sim7
jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim7/res_1_1 sim7_res_1_1
jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim7/res_1_2 sim7_res_1_2
jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim7/res_1_3 sim7_res_1_3
jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim7/res_1_4 sim7_res_1_4
jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim7/res_1_5 sim7_res_1_5

jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim7/res_2_1 sim7_res_2_1
jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim7/res_2_2 sim7_res_2_2
jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim7/res_2_3 sim7_res_2_3
jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim7/res_2_4 sim7_res_2_4
jl1 collect_results.jl /scratch/midway2/mkolar/diffTest/sim7/res_2_5 sim7_res_2_5


#!/bin/bash

## collects the results of the coverage simulations from $SCRATCH/diffTest into
## sim*/results, where the summarize.jl scripts read them
##
## ./collect.sh [sim ...]
##
## collects sims 1-4, 6 and 7 unless given. JULIA (default: julia) and SCRATCH
## (default: $HOME/scratch) can be set in the environment.

JULIA=${JULIA:-julia}
SCRATCH=${SCRATCH:-$HOME/scratch}

cd "$(dirname "$0")"

sims=${@:-1 2 3 4 6 7}

for sim in $sims
do
case $sim in
1|2) numElem=3 ;;
*) numElem=5 ;;
esac
mkdir -p sim${sim}/results
for ip in 1 2
do
for iElem in $(seq 1 $numElem)
do
echo "sim${sim} ${ip} ${iElem} ..."
$JULIA --project=. collect_results.jl $SCRATCH/diffTest/sim${sim}/res_${ip}_${iElem} sim${sim}/results/sim${sim}_res_${ip}_${iElem}
done
done
done
35 changes: 18 additions & 17 deletions exper/collect_results.jl
Original file line number Diff line number Diff line change
@@ -1,26 +1,27 @@
# Collect the estimates of one setting of a coverage simulation (sims 1-4, 6 and
# 7) into a matrix with one column per repetition:
#
# julia --project=. collect_results.jl resPath outPrefix
#
# reads resPath_rep.jld2 for rep = 1, ..., NUM_REP (default 1000) and writes
# outPrefix.jld2. Missing repetitions are skipped with a warning.
using DiffPrecTest
using JLD
using JLD2

resPath = ARGS[1]
outPrefix = ARGS[2]

# resPath = /scratch/mkolar/diffTest/sim1/res_2_2_$(rep).jld
NUM_REP = parse(Int, get(ENV, "NUM_REP", "1000"))

NUM_REP = 1000
results = Array{Any}(undef, 9, NUM_REP)
estimates = Vector{Any}[]
for rep=1:NUM_REP
global results
fname = "$(resPath)_$(rep).jld"

# try
file = jldopen(fname, "r")
res = read(file, "est")
close(file)

results[:, rep] = deepcopy(res)
# catch
# @show fname
# end
fname = "$(resPath)_$(rep).jld2"
if isfile(fname)
push!(estimates, load(fname, "est"))
else
@warn "missing $fname"
end
end
results = reduce(hcat, estimates)

@save "$(outPrefix).jld" results
@save "$(outPrefix).jld2" results
25 changes: 12 additions & 13 deletions exper/qq_plot.jl
Original file line number Diff line number Diff line change
@@ -1,23 +1,22 @@
using Distributions, StatsBase
using PyPlot
using Distributions
using Plots


function qqplot(z; color="black")
# normal Q-Q plot of the sample z added to plt
function qqplot!(plt, z; color=:black)
n = length(z)

grid = [(1 / (n + 1)):(1 / (n + 1)):(1.0 - (1 / (n + 1)));]
grid = (1:n) ./ (n + 1)

qz = quantile(z, grid)
qd = quantile.(Ref(Distributions.Normal()), grid)
qd = quantile.(Normal(), grid)

lims = 3.290
x = range(-lims, stop=lims)

plot(x, x, color="grey", linestyle=":", linewidth=.25)
scatter(qz, qd, s=.75, color=color)
plot!(plt, [-lims, lims], [-lims, lims]; color=:grey, linestyle=:dot, linewidth=0.5, label="")
scatter!(plt, qz, qd; markersize=1, markerstrokewidth=0, color=color, label="")

xlim([-lims, lims])
ylim([-lims, lims])
xlims!(plt, -lims, lims)
ylims!(plt, -lims, lims)

nothing
end
plt
end
Loading
Loading