nignatiadis / Empirikos.jl

Empirical Bayes estimation and inference in Julia.
MIT License
17 stars 0 forks source link

Support for Gaussian scale mixture class with a pointmass at 0 #11

Open DongyueXie opened 2 years ago

DongyueXie commented 2 years ago

Hi Nikos,

Thanks for providing a software implementation of the nonparametric EB confidence interval method.

I tried to run the simulation in simulation_ebci.jl with spiky, lfsr setting. I modified the Gaussian scale mixture class, gcal_scalemix, to include a pointmass at 0, as

gcal_scalemix = GaussianScaleMixtureClass([0, 0.1000000, 0.1100000 , 0.1210000, 0.1331000 ,... ,11.7390853, 12.9129938, 14.2042932, 15.6247225,])

However the confidence intervals of lfsr obtained do not seem quite right as shown below:

image

The shaded orange area is the confidence interval when using SN without point mass, and the dashed lines are upper and lower bounds of confidence intervals using SN with point mass at 0. The upper bounds are pretty close while the latter lower bound gets even lower. But the upper bound is expected to increase after adding the pointmass at 0.

Similarly, when using SN with point mass at 0, I tried targets = Empirikos.PosteriorProbability.(StandardNormalSample.(ts), Interval(0, 0)) and hope to get confidence intervals for local false discovery rate(lfdr). But the returned lower, upper bound are all 0's for all target z values.

So I was wondering if the software implementation supports Gaussian scale mixture class with a pointmass at 0? Or if there's any way to obtain the confidence interval when adding the pointmass?

Thank you! Looking forward to your reply.

nignatiadis commented 2 years ago

Dear Dongyue,

Thank you for reaching out about this! I had actually been worried/thinking about this just a bit before you opened this issue.

In our paper, we do not consider the case of a point-mass at zero (and indeed, to me the beauty of the local false sign rate is that it makes sense without requiring exactly null effects). The consequence however is that the current code does not correctly account for point masses (it should throw an error on Normal(0,0)).

I will try to add proper support for point masses in the code, however, in the meantime the following "hack" should work: Use Interval(-eps(), nothing) for the Local False Sign Rate and Interval(-eps(), +eps()) for the Local False Discovery Rate, i.e.,

targets_lfdr = Empirikos.PosteriorProbability.(StandardNormalSample.(ts), Interval(-eps(), +eps()))
targets_lfsr = Empirikos.PosteriorProbability.(StandardNormalSample.(ts), Interval(-eps(), nothing))

Let me know if this works OK for now and also if you have any questions.

Best, Nikos

DongyueXie commented 2 years ago

Hi Nikos,

Thanks for your reply.

Your suggested method of using eps() seems to work when combining with GaussianScaleMixtureClass with a point-mass, i.e. set gcal_scalemix = GaussianScaleMixtureClass([0, 0.1000000, 0.1100000 , 0.1210000, 0.1331000 ,... ,11.7390853, 12.9129938, 14.2042932, 15.6247225,])

Then the confidence interval of lfsr looks more reasonable:

image

If not adding that 0 in gcal_scalemix, then using eps() does not change the results much.

Also the current code does not return error when using Normal(0,0), either for generating effects or in above gcal_scalemix with a 0 sd.

nignatiadis commented 2 years ago

Hi Dongyue,

That sounds right to me! So to recap:

Best, Nikos