negative values in confidence interval of multi-taper coherence

[Xunius]

First of all, I still need to read more carefully the references, so I might be wrong.

In the Multitaper coherence estimation tutorial, the confidence interval are computed (t975_limit and t025_limit), but they are not printed or visualized in any way later. It turns out that t025_limit contains many negative values, but coherence is constrained to be within [0, 1].

Is there anything going wrong?

[miketrumpis]

I thought I'd weigh in here (hi @arokem !).. the problem is with computing t-distributed confidence intervals, scaled by the jackknife standard error:

t025_limit = coh_mat_xform + dist.t.ppf(.025, K - 1) * np.sqrt(coh_var)
t975_limit = coh_mat_xform + dist.t.ppf(.975, K - 1) * np.sqrt(coh_var)

But, as you correctly point out, coherence is [0, 1] (very much not Normal). I think the typical thing to do here is use the Fisher transform (tanh & arctanh). That would be something like

upper_ci = tanh(arctanh(mu) + t_ci * SE)
lower_ci = tanh(arctanh(mu) - t_ci * SE)

arctanh makes the [0, 1] value Normal-ish, so you can use t-distributed CIs. Then tanh brings you back to the unit scale. I believe something similar is done for PSD CIs (operating on the log-PSD as a Normal-ish quantity).

edit: I guess you'd want to compute the Jackknife SE with arctanh values as well (??)

[miketrumpis]

On second read.. now that I'm looking closer, the script is already doing the normalize / un-normalize steps. @Xunius are you saying the upper/lower CIs are still outside the [0, 1] range? This block literally outputs tanh values

utils.normal_coherence_to_unit(t025_limit, 2 * K - 2, t025_limit)
utils.normal_coherence_to_unit(t975_limit, 2 * K - 2, t975_limit)

[Xunius]

@miketrumpis Yeah but tanh() is within [-1, 1]

[miketrumpis]

edit: disregard, bad idea.. I'll give this some more time

Good point! Maybe offset the basis?

def coh_to_norm(x):
    return np.arctanh((x + 1) / 2)

def norm_to_coh(y):
    return (np.tanh(y) + 1) / 2

[miketrumpis]

Those transforms weren't quite right. This appears to handle the tails appropriately (see gist)

Also, these are simplified (not scaled by DoF)

def coh_to_norm(x):
    return np.arctanh(2 * x - 1)

def norm_to_coh(x):
    return (np.tanh(x) + 1) / 2

[Xunius]

The new norm_to_coh() and coh_to_norm() functions look reasonable.

I'm going through the reference Thompson, DJ (2007) Jackknifing multitaper spectrum estimates. IEEE Signal Processing Magazing. 24: 20-30. I think the current code does a faithful translation of the Equation (24) in the paper. But I still have a couple of questions:

1. The paper mentioned "see [4] for details references". But I couldn't find a copy of the reference [4]: "Jackknife error estimates for spectra, coherences, and transfer functions". This is the only thing I can find on google. Does anyone find a copy of this?
2. The Thompson 2007 paper mentioned magnitude-squared-coherence (MSC), but Equation (24) gives the a definition of coherence that is different from the nitime definition in that Equation (24) coherence is the sqrt of MSC. So I wonder whether we are using the wrong equation? Does this affect the scaling factor 2*K-2 used in normalize_coherence() and normal_coherence_to_unit()?

[miketrumpis]

I'm seeing eq 24 as complex, but I agree: it does look like the 2007 sig proc magazine paper is suggesting to use the arctanh transform of the modulus of the (complex) coherency find a jackknife standard error.. hard to say without that book (it's not even listed at the Duke library). Anyway, that inference is consistent with how jackknifed_coh_variance works.

I believe MSC is more common in engineering, and the not-squared magnitude coherence is more common in neurophys. I think the scaling factor is appropriate for the not-squared.

But the question is probably that the "one-sided" Fisher transforms are a total hack, and are they valid? Probably not! I hope they're better conditioned and numerically stable, though. Not a statistician, but my hunch is distribution symmetry is probably more important than normality when using resampling methods.

There is a body of literature on this, so I hope there's a real answer. I am a little band-limited at the moment. If you'd like to dig in more, there should be other sources to find those parametric estimator stats related to arctanh (my memory is that's a "classic" stats finding). This might be helpful

Bokil J neuro meth 2007

ps: I'm scanning this and assuming that there is an absolute value between eq 3 (definitely complex) and eq 4, because they distinguish coherence from coherency in their language.

[Xunius]

The Equation (24) in Thompson 2007 paper is the not-squared coherence, and it is inconsistent with the term MSC. So it is either the MSC term is used incorrectly or the equation is wrong.

The context where these jackknife confidence intervals are computed, i.e. the MTCoherenceAnalyzer, uses a definition of coherence that is magnitude squared. So whatever convention you neurophys people prefer, this is inconsistent.

Regarding to the transformation, it may be helpful to if I could have a read of the citation Thompson referred to (Jackknife error estimates for spectra, coherences, and transfer functions). Do you happened to have a copy of that?

I found another package, and it seems to give quite different results in some estimates (MT coherence, spectrum) than nitime. Not sure exactly why though.