search

NeuroStat / neuRosim

An R package to simulate fMRI Data Including Activated Data, Noise Data and Resting State Data
GNU General Public License v2.0
11 stars 6 forks source link

physnoise.r does not satisfy SNR when SNR > 1 #5

Closed JuliaBeitner closed 5 years ago

JuliaBeitner commented 5 years ago

Dear maintainers, I have a question regarding the function creating physiological noise (i.e., physnoise.R). In line 15, the additional noise is drawn from a normal distribution with a mean of 0 and SD of 1, while the physiological noise itself is rescaled to fit the given parameter sigma. noise <- array(rnorm(prod(dim)*nscan, 0, 1), dim=c(dim, nscan)) + HRweight*HRdrift In cases where the given SNR is larger than 1, sigma becomes smaller than 1 (e.g., SNR = 10, sigma = 0.1, if average signal is 1). The output, however, will then not have a SD close to sigma but to 1 because of the additional noise from the normal distribution. Also, in this case the additional noise is stronger than the physiological noise itself. When using the simTSfmri.R function with mixture noise (or also the simTSrestingstate.R and simVOLfmri.R), the given SNR does not reflect the SNR of the simulated data anymore, if the SNR exceeds 1.

The following solutions came to my mind: 1) Remove the additional noise part and only keep the physiological noise. 2) Rescale the original output with noise <- noise * sigma/sd(noise) (line 32 from tasknoise.R). 3) Change the SD of the normal distribution in such a way, that the proportions of additional and physiological noise are kept constant and the final output satisfies the prespecified sigma.

In all three cases, the function would now satisfy given SNRs, but they all model physiological noise in a different way.

What would you advise when using a mixture of all noise types with the default weights (30% system (white) noise, 30% temporal noise, 1% low frequency drift, 9% physiological noise, and 30% task noise)?

jdgryse commented 5 years ago

Dear Julia,

Thank you for pointing out this issue. We made a new branch in the neuRosim repository, SNR, in which we updated the code of several files, so that the SNR of the simulated data now reflects the value that is used as input in the simVOLfmri() function for example. This will also be the case when adding physiological noise. In the script physnoise.R, Gaussian noise is still added to the purely physiological noise, however, an issue with the code for how the latter was constructed is now resolved.

Would you be willing to verify whether the obtained SNR is indeed equal to the SNR that you asked for in your code? If so, we can then merge the SNR branch with the master branch. Thank you!

JuliaBeitner commented 5 years ago

Dear maintainers,

Thank you very much for your reply. I checked, and the SNR is now correct and equal to the SNR specified. So the issue can be closed. However, there was a new bug in tasknoise.r (please see https://github.com/NeuroStat/neuRosim/pull/6). Besides that, everything was working fine and as it was supposed to be. Thank you again!