issue with knee model

[ZoeD2A]

Hi,

I have a problem with the fm.fit. All my knee values are very high. I have checked my power spectrum and it looks OK (no artefacts), but as soon as I run the Fooof analysis and fit, the knee values are totally inconsistent.
I tried with the FOOOF and then with the SpectralModel and it gives the same thing. Can anyone help me?

Here is a part of my code:

Parameters

sampling_rate = 250
min_peak_height = 0.15
max_n_peaks = 6
min_freq = 4
max_freq = 75
peak_width_limits=[2,25]
frequency_range = (min_freq, max_freq)

for index, row in filtered_df.iterrows():
if isinstance(row["Ch 1 data"], (list, np.ndarray)) and len(row["Ch 1 data"]) > 0:
freqs, powers = {}, {}

                # Compute Welch power spectrum for all channels        
                for ch in range(1, 5):
                    channel_data = row[f"Ch {ch} data"]
                    freqs[ch], powers[ch] = compute_spectrum(channel_data, fs=sampling_rate, method='welch',avg_type='mean', nperseg=500, scaling='density')    

                    ## added to limit the frequency range
                    freq_mask = (freqs[ch] >= min_freq) & (freqs[ch] <= max_freq)
                    
                    if len(freqs[ch]) > 0 and len(powers[ch]) > 0:
                        filtered_df.at[index, f"Ch {ch} frequencies"] = freqs[ch][freq_mask] 
                        filtered_df.at[index, f"Ch {ch} power_spectrum"] = powers[ch][freq_mask] 
                        fm = SpectralModel(peak_width_limits,aperiodic_mode='knee', min_peak_height=min_peak_height, max_n_peaks=max_n_peaks, verbose=False)

                    try:
                        print(f"Running FOOOF on Channel {ch}, Row {index}")
                        fm.fit(freqs[ch][freq_mask], powers[ch][freq_mask], frequency_range) 
                        fm.report() 

FOOOF results if fm.has_model:

                            print(f"Model fit successful for Channel {ch}, Row {index}")
                            aperiodic_params = fm.get_params('aperiodic')
                            filtered_df.at[index, f"Ch {ch} offset"] = fm.get_params('aperiodic')[0]
                            if len(fm.aperiodic_params_) > 2:
                                filtered_df.at[index, f"Ch {ch} knee"] = fm.get_params('aperiodic')[1]
                                filtered_df.at[index, f"Ch {ch} log_knee"] = math.log10(max(aperiodic_params[1], 1e-6))
                            filtered_df.at[index, f"Ch {ch} exponent"] = fm.get_params('aperiodic')[2]
                            filtered_df.at[index, f"Ch {ch} r_squared"] = fm.get_params('r_squared')
                            filtered_df.at[index, f"Ch {ch} error"] = fm.get_params('error')
                            filtered_df.at[index, f"Ch {ch} peaks"] = fm.get_params('peak_params')[0]
                            print(f"Results for Channel {ch}, Row {index}: Aperiodic parameters: {aperiodic_params}")

And the output for FOOOF analysis:

Thank you for your help!

[voytek]

Can you share some of the spectra? It helps to see what they look like.

[ZoeD2A]

Hi @voytek,
Thank you for you answer.
Here some examples of my spectra.

[voytek]

Hmmm, can you also share the resulting image with the spectrum and the model that you get from fm.report()?

[ZoeD2A]

The direct results from fm.report show for example:

But when I apply FOOOF analysis separately to one of the signals, I have for example:

[voytek]

Do you get the wildly large knee values in both cases? The fitting range you're using is very different between the two sets and I'm wondering if that's somehow causing this.

[ZoeD2A]

Sorry, I don't really understand your answer. But here are two more examples:

My code to plot the signal in my first channel:

file_path = os.path.join(out_dir, 'patient_power-spectrum_FOOOF-TEST.csv')
data = pd.read_csv(file_path)

index = 3
day_night_value = data['Day/Night'].iloc[index]
print(f"Day/Night for {index}: {day_night_value}")

def parse_list(string):
return np.array([float(i) for i in ast.literal_eval(string)])
frequencies_1 = parse_list(data['Ch 1 frequencies'].iloc[index])
power_spectrum_1 = parse_list(data['Ch 1 power_spectrum'].iloc[index])

freq_mask = (frequencies_1 >= 4) & (frequencies_1 <= 75)
frequencies_filtered = frequencies_1[freq_mask]
power_spectrum_filtered = power_spectrum_1[freq_mask]
log_power_spectrum_1 = np.log10(power_spectrum_1 + 1e-10)
log_power_spectrum_filtered = log_power_spectrum_1[freq_mask]

frequency_range = (4, 75)
fm = SpectralModel(aperiodic_mode='knee', min_peak_height=0.15, max_n_peaks=6, verbose=False)
fm.fit(frequencies_filtered, power_spectrum_filtered)
fm.report()

plt.figure(figsize=(12, 8))
plt.plot(frequencies_filtered, log_power_spectrum_filtered, label="Log Raw Spectrum (Ch 1)", color='orange', alpha=0.6)
fm.plot(ax=plt.gca(), plot_peaks='shade', plot_kwargs={'color': 'blue', 'linestyle': '--'})
#plt.plot(frequencies_1, power_spectrum_1, log_power_spectrum_1, label="Log Power Spectrum (Ch 1)", color='blue', alpha=0.6)
#fm.plot(ax=plt.gca(), plot_peaks='shade', plot_kwargs={'color': 'red', 'linestyle': '--'})

for peak_params in fm.get_params('peak_params'):
peak_freq, peak_power, _ = peak_params
log_peak_power = np.log10(peak_power + 1e-10)
plt.annotate(f'{peak_freq:.1f} Hz', (peak_freq, log_peak_power), textcoords="offset points", xytext=(0, 10), ha='center', color='darkred')

plt.title(f"Power Spectrum and FOOOF Fit (Index {index}, Day/Night: {day_night_value})")
plt.xlabel("Frequency (Hz)")
plt.ylabel("Log(Power)")
plt.legend(["Log Raw Spectrum (Ch 1)", "FOOOF Fit", "Peaks"])
plt.grid(True)
plt.show()

And two examples of my results:

the power spectrum:

                                POWER SPECTRUM MODEL                                       
                                                                                              
                    The model was run on the frequency range 4 - 75 Hz                        
                             Frequency Resolution is 0.50 Hz                                  
                                                                                              
                     Aperiodic Parameters (offset, knee, exponent):                           
                                 3.7060, 140.6524, 2.3725                                     
                                                                                              
                                   4 peaks were found:                                        
                            CF:  23.27, PW:  0.090, BW:  8.62                                 
                            CF:  40.83, PW:  0.159, BW:  0.50                                 
                            CF:  45.72, PW:  0.109, BW: 12.00                                 
                            CF:  50.31, PW:  0.151, BW:  0.50                                 
                                                                                              
                                 Goodness of fit metrics:                                     
                                R^2 of model fit is 0.9909                                    
                                Error of the fit is 0.0481 

The power spectrum:

POWER SPECTRUM MODEL

                    The model was run on the frequency range 4 - 75 Hz                        
                             Frequency Resolution is 0.50 Hz                                  
                                                                                              
                     Aperiodic Parameters (offset, knee, exponent):                           
                                 2.9318, -11.2267, 2.0876                                     
                                                                                              
                                   4 peaks were found:                                        
                            CF:  28.08, PW:  0.147, BW:  9.30                                 
                            CF:  36.02, PW:  0.219, BW:  2.24                                 
                            CF:  42.33, PW:  0.251, BW:  0.66                                 
                            CF:  44.92, PW:  0.128, BW: 12.00                                 
                                                                                              
                                 Goodness of fit metrics:                                     
                                R^2 of model fit is 0.9948                                    
                                Error of the fit is 0.0369 

Thank you for your help,

[voytek]

Right now you're fitting from 4-75 Hz. What if you change it to 1-60 Hz?

[ZoeD2A]

It's better but I always have strange values. I don't understand why I have a high knee value.
For example,

For the channel 1:
POWER SPECTRUM MODEL

                    The model was run on the frequency range 1 - 60 Hz                        
                             Frequency Resolution is 0.50 Hz                                  
                                                                                              
                     Aperiodic Parameters (offset, knee, exponent):                           
                                 3.6490, 606.8577, 2.3896                                     
                                                                                              
                                   5 peaks were found:                                        
                            CF:   2.96, PW:  0.703, BW:  1.49                                 
                            CF:   5.43, PW:  0.310, BW:  2.85                                 
                            CF:  19.24, PW:  0.228, BW:  5.19                                 
                            CF:  29.34, PW:  0.154, BW:  2.49                                 
                            CF:  32.64, PW:  0.218, BW:  1.70                                 
                                                                                              
                                 Goodness of fit metrics:                                     
                                R^2 of model fit is 0.9928                                    
                                Error of the fit is 0.0362 

[voytek]

@ryanhammonds do you have any recommendations or suggestions here?

[ryanhammonds]

This is related to #224 and #235. The knee parameter as it is right now is not very interpretable. To convert it to a frequency, use:

knee_freq = knee**(1./exponent)

For the params reported above it gives 14.6 Hz knee freq, e.g.606.8577**(1/2.3896). My fork here directly optimizes for knee_freq.

[ZoeD2A]

Thank you for your answers!

[TomDonoghue]

As far as I can tell, this is a question that has been addressed - so I'm going to close this now - but feel free to reopen / follow up if there's more to this to check in on!