Introduction

This demo reproduces the ordered hypothesis testing simulation in the paper:

The demo was run using R version 3.2.0.

The demo considers an ordered hypothesis testing problem with \(n\) hypotheses, of which \(n_1\) contain true signals. We produce an ordered list of p-values as follows:

We then run accumulation tests, specifically, we text the HingeExp method (with parameter \(C=2\)), the SeqStep and SeqStep+ (Barber and Candès 2014) (each with parameters \(C=2\)), and the ForwardStop method (G’Sell et al 2013).

Setup

We first define some functions for running the simulation and for plotting the results. (Code for this setup hidden in the output, but is visible in .Rmd file)

Setup_simulation_functions()
## [1] "Simulation functions defined."
Setup_simulation_plotting_functions()
## [1] "Simulation plotting functions defined."
source('accumulation_test_functions.R')

Parameters for the simulations:

hfuns=c(create_SeqStep_function(C=2),create_SeqStep_function(C=2),create_HingeExp_function(C=2),create_ForwardStop_function())
numerator_pluses=c(2,0,0,0)
denominator_pluses=c(1,0,0,0)
names=c('SeqStep+ (C=2)','SeqStep (C=2)','HingeExp (C=2)','ForwardStop')
mu_low=2;mu_high=3; 
n=1000; ntrials=100; 
pr_list=c(0.2,0.1) # proportion of true signals
alphaseq=(2:10)/40
seeds=1:ntrials

Simulations part 1

Settings: \(n=1000\) hypotheses, \(k^\star=200\) true signals

First, we run each method with target FDR level \(\alpha\) (ranging in \(\alpha=\{0.05,0.075,0.1,\dots,0.25\}\)). For each method and each \(\alpha\), we record the actual FDR attained, and the power to detect signals, averaged over 100 trials.

pr=pr_list[1]
#pdf(paste0('FDR_vs_power_plot_n_',n,'_kstar_',n*pr,'.pdf'),14,14)
par(mfrow=c(2,2))
for(i in 1:2){for(j in 1:2){
  mu1=c(mu_low,mu_high)[i];mu2=c(mu_low,mu_high)[j]
  temp=CompEffectFunc(n,pr,mu1,mu2,hfuns,numerator_pluses,denominator_pluses,seeds, alphaseq)
FDR_and_power_plot(temp[[2]],temp[[1]],names,cols=c('black','blue','purple','red'),pchs=15:18,alpha=alphaseq, alpha_display_div=2, title=get_titles(mu_low,mu_high)[i,j],show_legend=(i==1 && j==1))
}};#dev.off()