Getting voxelwise maps for logL (and other attributes of LikelihoodModelResults)

[jdkent]

To get residual maps nilearn/nistats#410, we added a method called _get_voxelwise_model_attribute, which has broader use than returning residuals, such as "potentially" returning log-likelihood and other attributes of RegressionResult and LikelihoodModelResults.

However, there was some difficulty getting a voxelwise map of the logL results, with discussion/investigation outlined here and potential solution outlined here

Let me know if the solution looks reasonable, and I can open a pull request.
(and perhaps outline which attributes we would want to make explicit methods for in addition to residuals, predicted, and r_square in the FirstLevelModel class. For example, including:

• logL
• SSE
• MSE
• ???

[adelavega]

I'd also like to see a degrees of freedom map

[jdkent]

Do degrees of freedom vary by voxel? My understanding is that every voxel should have the same number of time points and every voxel should have the same model applied to it. (here is the calculation) I could be looking in the wrong place...

Would you like an single value attribute representing degrees of freedom be available from the FirstLevelModel class if degrees of freedom do not vary between voxels?

[bthirion]

Indeed, degrees of freedom are constant across voxels.

[bthirion]

@jdkent
Yes, please open a PR.

[adelavega]

Great, in some other software it can vary by voxel if its estimated, so no need for that, then.

[tyarkoni]

@bthirion, is there a principled reason to store voxels under separate AR1 models, as opposed to the other way around? It seems more intuitive to me to store all the attributes (R^2, logL, etc.) as 1d arrays of length n_voxels than to have variable length arrays under each of the AR1 bins. To ensure the mapping is preserved, one could add an ar_model attribute that provides each voxel's index into the list of ARModels. That seems like a cleaner API and would eliminate the _get_voxelwise_model_attribute step. But I may be missing something about how the computations are performed internally that would make this approach less practical.

[bthirion]

Thx for you suggestion

• First by default I really don't want to store R^2 logl and all these quantities because in all use case I know, people don't need them, and they consume memory, which ultimately an issue. We compute and store them if needed but not as a default behavior.
• From what I understand, you would like to store all the statistics in arrays of size (feature, dimension), which is indeed clearer to introspect the object.
  This is however a bit clumsy when computing contrasts, because we still rely on the underlying regression objects.
  See e.g.
  https://github.com/nilearn/nistats/blob/master/nistats/contrasts.py l.73
  So, I think that this may be changed (and probably should), but this requires some refactoring.