-
Notifications
You must be signed in to change notification settings - Fork 64
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
[ENH] - Add autocorrelation related functionality #331
Conversation
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Looks good! I left one comment about some confusion of units of tau (seconds vs samples).
Future ideas not totally relevant to this PR:
- If this is PR is in scope, a minimal specparam version (e.g. only the simple_ap_fit method) may also be.
- It may be possible to find equivalent, 1:1 forms, here is the spectral form that is equivalent to acf exponential decay. In these cases, specparam could contain a "parameters_to_acf" like method to handle conversion.
- Given equivalent forms, when is it more accurate to optimize parameters in ACF vs PSD space?
Results of fitting the function to the autocorrelation. | ||
""" | ||
|
||
return scale * (np.exp(-timepoints / tau) + offset) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This tau is in samples. Could make a note in docstring that tau isn't in seconds (what I expected). Or the func could accept fs as a separate argument, then do the multiplication in the func. Here is an example:
import numpy as np
from neurodsp.aperiodic.autocorr import *
fs = 1000
tau = 0.015
lags = np.linspace(0, 100, 1000)
corrs = exp_decay_func(lags, tau * fs, 1, 0) # must scale tau by fs to get correct results
params = fit_autocorr(lags, corrs, fs)
print((params[0], tau)) # returns (0.01499988686748319, 0.015)
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
You're right - there's an inconsistency here.
fs
isn't a required inpit to fit_autocorr
, so if not passed, then it's consistent (this make sense, I think):
corrs = exp_decay_func(lags, tau, 1, 0)
params = fit_autocorr(lags, corrs)
print((params[0], tau)) # returns (0.014996887598931195, 0.015)
But if you do pass fs
into fit_autocorr
, then the timepoints get rescaled, but the tau value does not, leading to the issue.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Another way to write what you point out is that if fs
is passed into the fit function, then the timepoints
should be in time not in samples, and so recomputing the values consistently should look like (this is equivalent to multiplying the tau):
corrs = exp_decay_func(lags / fs, tau, 1, 0)
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Okay, the core of this issue is basically it's not a good idea to optionally take fs
as an input in fit_autocorr
and if so divide it out, because then it makes the values of the fit parameters different from the values of the inputs. I'll drop this, and then one can pass in timepoints as either samples or seconds, and it will work either way (and be consistent).
Thanks for the review @ryanhammonds - I've pushed some updates that I think cleans things up and addresses the issue you addressed. I'm going to merge this in now, since I'm trying to consolidate some branches to use the new version in a project - if there is anything to follow up on here we can do another PR and make sure to address it before tagging the next version! I also agree that once the next specparam version comes together with the additiona models, we can revisit the alignment between methods and see what should be added where! |
When doing the Aperiodic Methods project, I ended up using some autocorrelation related funtionality - notably measuring AC decay times, and fitting functions to AC curves - that I think might as well live in neurodsp.
I adapted the functions here from @rdgao's field echos / knee project