# Should I use an average ECDF?

This relates to a previous question of mine which didn't gain many responses, perhaps because it wasn't very clear nd well written. I hope this time I will be more accurate and get your much appreciated assistance.

I am analyzing results of a biological experiment. The results given as a single value ( non-negative integer) per genomic position. I am interested in valleys, or local minimam over this series of values.

I wish to control the false positives rate and get the significance for each local minima. I can shuffle the raw data which was used to produce the data.

So what I do is to shuffle the raw data, create the new series of values,search for all local minima and keep their values.

Now, I have something like this:

data_set  local_minima_values
=============================
true_data 4 9 1 27 12 0 0 2 5 32 0 1 5 70 2
sim_1 14 25 94 59 32
sim_2 52 0 14 74 82 12 54
...


Note the number of local minima naturally varies between simulations.

So, my idea was to calculate an ECDF for each simulation and then combine those ECDFs into a single "average ECDF" which represents the null hypothesis. Then, I can assign a p-value for each local minima from the true data, and see how significant ('surprising') it is.

My questions are:

1. Does this make sense?
2. How do I create an average ECDF? I can't just merge the values from all simulation together and get and ECDF for this merged set, since the number of minima found in each simulation differs, and I think all simulations should have the same contribution to the average ECDF, or am I wrong?
3. How should I take the number of simulations (shuffles) into account?

Thanks,

Dave

p.s. I'm working with R.

To average the ECDFs, I'd do something like:

impute_resolution = 1e3
values_to_impute = seq(
min(my_data$true_data) , max(my_data$true_data)
, length.out = impute_resoluton
)

ecdfs = matrix(NA,nrow=length(values_to_impute))

for(i in 1:(ncol(my_data)-1)){ #assumes column 1 is true_data
this_ecdf = ecdf(my_data[,i+1])
ecdfs[i,] = this_ecdf(values_to_impute)
}

mean_ecdf = colMeans(ecdfs)
plot(
x = values_to_impute
, y = mean_ecdf
, type = 'l'
)


David, before discussing the implementation, I'd discuss question #1. I think the approach "makes sense", if you mean by this that the approach is intuitive. However, from the perspective of bootstrapping, the approach probably won't work. One well-known failure of bootstrap is in fact in the case of maxima and minima. The intuition behind this is that the estimators for these statistics have very large variance: a nonparametric approach based on the empirical distribution doesn't help estimate these statistics. Most approaches are parametric. If interested I can suggest references on the latter.

A reference on the failure of boostrap for extremal statistics is Efron and Tibshirani's very own "An Introduction to the Bootstrap". Check page 81.