Morlet wavelet transform normalization

[Moo-Marc]

Hello!

While working on this issue and the TF tutorial following my previous update to spectrum normalization, I found that the wavelet transform itself introduces a 1/f power factor, meaning that if for example you simulate 3 sine waves at different frequencies with equal amplitudes, our wavelet transform implementation (with the option 'none' for normalization) will produce 1/f power. See the first line in the figure below, where the second column is the average "spectrum" from 2 to 6 seconds computed by averaging the TF in time, similar to when we right click on the TF and select "spectrum". The third column shows the actual Welch spectrum for this simulated data.

The second line shows the same data with my "corrected" 1/f compensation to the power spectrum, while the 3rd row shows Matlab's Morlet continuous wavelet transform (cwt), without further processing (no explicit 1/f compensation).

The issue is not that this additional 1/f scaling of our implementation is "wrong", but rather that it was not really documented or explained anywhere, and that previously, the 1/f^2 compensation was actually compensating for both this transform effect, and the "real" neural 1/f spectrum. So for the wavelet transform, this factor could be considered ok, if confusing, but for the others (Hilbert and PSD), I think my correction from 1/f^2 to 1/f was appropriate.

Now, I imagine that the logic behind these wavelet coefficients dropping in power is that the power actually gets spread over an increasingly wide frequency range (wider peaks), and probably the total power is conserved. However, for viewing a TF decomposition of the data, I would expect scaling the coefficients to maintain peak power (instead of total power) is more appropriate.

I would therefore suggest renaming the options for the wavelet transform:
Spectrum normalization:

1. Conserve total power (previous 'none' option - peaks have an additional 1/f scaling, in grey)
2. Conserve peak power (didn't exist before my previous change)
3. 1/f compensation (previous '1/f' option which applied 1/f^2)

This way, "1/f compensation" will consistently refer to the neural effect throughout Brainstorm.

Thoughts?

[sbaillet]

Great catch, Marc -- thanks !
I like your idea but I am unclear how the "Conserve Peak Power" scaling would work? Indeed, how would you determine the actual peak power. I suppose this was explained in your previous change, as you mention above, so sorry if missed it but please provide a pointer.

[Moo-Marc]

how would you determine the actual peak power

I just mean it compensates for the 1/f peak amplitude scaling that our wavelet transform introduces: 2nd row in my figure. It doesn't compensate for neural 1/f.

[Moo-Marc]

From @dimitriospantazis :

I generally have not been using spectral normalization in my work; instead I have been using the common baseline normalization. But here is my take.

You are right that there are two 1/f normalizations: one for the neural 1/f effect (which applies to all power transformations), and one 1/f specifically for the Morlet wavelet scaling.

Also, something to note is that the continuous Morlet wavelet is not an orthonormal transformation so the power from the time domain to the spectral domain is not preserved. This means the Frobenius norm of the time domain and the spectral domain will be different. But I think the relative power of two sinusoidal signals will remain equal in the wavelet domain. Or in other words, the two Gaussian kernels in your plot below (centered in 20Hz and 50Hz) will have equal integrals.

In general we care about relative changes of power within the same frequency band. In fact, I am not sure we ever need to explicitly compare the power between frequency bands (e.g. say that alpha power is stronger that gamma power). Thus the 1/f is largely for visualization. Please correctly if I am wrong.

Let me make some observations for the Brainstorm implementation:

(1) Do we need the ‘spectrum normalization’ for reasons other than visualization? If not, then it appears more appropriate to include it in the visualization controls (next to ‘Power’, ‘Magnitude’, ‘Log(power)’) rather than as a separate process.

(2) I have noticed that the standardization process provides two options: ‘spectrum normalization’, which only changes a flag in the file, and ‘baseline normalization’, which is permanent. Maybe this is a bit counterintuitive.

(3) If the ‘spectrum transformation’ is not permanently applied to the data, then I assume it is because we want to keep it reversible. But Brainstorm does not allow me to change it once applied. See (*) below.

(4) For the Morlet wavelet transformation, the process window already has a 1/f transformation called called ‘Spectral flattening: Multiply output power values by frequency’. If we introduce a similar transformation in the ‘spectrum normalization’, we need to use consistent names.

(5) I am not entirely sure a (1/f)^2 transformation preserves the peak power. This is not obvious in the figure above where one peak is 1.1 and other is 1. Can you please verify? Even if that is the case, two sinusoidal signals with closer frequencies (or three sinusoidal signals, etc) will yield variable peaks.

(6) I think the names ‘conserve total power’, ‘conserve peak power’, and ‘1/f compensation’ are not helpful. I am afraid they would cause confusion because (a) these names will need to differ depending on whether the data came from the Morlet wavelet or another transformation, (b) at that stage of analysis, the input data is power data, and you do not conserve their properties (the conservation is with respect to the original time series, but feels rather disconnected at that stage), and (c) the conservation is a rather abstract idea, between sinusoidal signals across frequencies. I would rather have ‘1/f compensation’ and ‘1/f^2 compensation’ as it is now. In general, researchers and especially neuroscientists are familiar with the 1/f transformation.

If I were to propose some suggestions, perhaps these steps:

-Simple approach: Keep names ‘1/f compensation’ and ‘1/f^2 compensation’. Also, for consistency, rename the Morlet process to read ‘1/f compensation (multiply power by frequency)’ instead of ‘Spectral flattening: Multiply output power values by frequency’

-Complicated approach: If only for visualization, remove processes (‘spectral normalization’ from the process window and the Morlet process) and add Brainstorm controls for visualization only. If for processing too, make it a permanent transformation?

These are only some initial thoughts and a team feedback would be useful. For example, one could consider two check boxes with names ‘1/f compensation (brain)’ and ‘1/f compensation (Morlet wavelet)’ and allow the selection of one or both. Also, feel free to move my comments to the github links below is you think it is more helpful.

Best,

Dimitrios

[Moo-Marc]

Thanks Dimitrios,

Point 1 is a good idea, though I don't think we'd want to remove a feature that people might have been using. Because actually, "spectrum normalization" does change the data, it's not just a flag. That's why the process checks and won't perform it repeatedly. But we could have a new visualization feature that would allow changing it on the fly as well. We'd want to be careful that the two (applied to data vs visualization) is not confusing though. Maybe make it red like the visualization filter tab.

Regarding consistent naming, I agree. In my next PR I removed the name "spectral flattening" which I didn't like and it was called "spectrum normalization" in other places. But consistency is also why I wanted to keep "1/f compensation" to consistently refer to the brain feature only, and use a different name when referring to the Morlet wavelet transform. I don't have a great idea for naming, my initial suggestion was just based on what I observed. I suggest we think of something better. For now, I'll refer to it as the "wavelet scaling option", with either "decreasing scale" (our implementation, first line in my original figure) and "increasing scale" (second line in my figure). (The reason for this temporary naming will become clear with my new tests below.)

This new option would of course only apply to the wavelet and would only appear in the wavelet transform process. The tf normalization process offers the options '1/f compensation', and '1/f^2 compensation' for backwards compatibility.

I did another test with real brain data to see what would be the effect on peaks that aren't a pure sine wave. Since higher frequency peaks are wider, the overlap indeed makes them taller with the "increasing scale", and the background also increases. With the "decreasing scale", the background stays flatter, but the peaks decrease in the simulation.
Real data:

Simulation with equal amplitude sinusoids:

I computed the peak power ratio between the wavelet spectra and PSD spectra.
Real data, 10, 20, 30 Hz:

• "increasing scale" power / PSD power = 4.1239 6.5244 14.1439 (increasing)
• "decreasing scale" power / PSD power = 0.3924 0.3199 0.4646 ("flat")

Simulated sinusoids, 2, 20, 50 Hz:

• "increasing scale" power / PSD power = 3.0087 3.0999 3.1859 ("flat")
• "decreasing scale" power / PSD power = 1.5043 0.1550 0.0637 (decreasing)

Looking forward to further discussion, and better suggestions for naming the two wavelet scaling options.

[dimitriospantazis]

Hi Marc,

Thanks for looking into this and providing the new plots. At this stage I think a feedback from the rest of the team would be useful. How about we defer this discussion to our next Brainstorm zoom meeting.

Cheers,
Dimitrios

[ftadel]

Any update on this topic?

[ftadel]

I updated the tutorial to show these new changes (1/f of erroneous 1/f^2 previously).
The labels and results had changed sensibly, I fixed the text and most screen captures.
https://neuroimage.usc.edu/brainstorm/Tutorials/TimeFrequency

[Moo-Marc]

Thanks François. Sorry I didn't find time to work on this since we discussed it in a Wednesday meeting back in Feb. Sadly, I'm having difficulty finding what we had decided. @dimitriospantazis please correct me if I'm wrong, or we can quickly confirm at the next meeting:

1. Add visualization options for spectrum/TF normalization.
2. Remove spectrum normalization option from TF process, and probably also from spectrum process (but keep as hidden options same as they were before for script backwards compatibility).
3. In spectrum normalization process, which would now be less useful because of point 1, make 1/f^2 option not grey again, since it can be useful for visualizing TF plots.
4. Update tutorials!

[dimitriospantazis]

Marc, I think the summary above is accurate. Regarding point 3, we can keep the 1/f and 1/f^2 normalization as processes to simplify things, even though they are likely not useful any more since we only need them for visualization (point 1).

Best,
Dimitrios

[ftadel]

@Moo-Marc Updates on this topic?