Track and write out S0 and T2* estimate variance (and maybe covariance?) maps
Summary
scipy.optimize.curve_fit outputs covariance estimates. We could compile those into spatial maps indicating T2* estimate variance, S0 estimate variance, and T2*-S0 estimate covariance.
I was just looking into something similar. We are currently saving the RSME of the S0 and T2* estimates as a quality metric. A potentially more useful measure is R^2 (i.e. the variance explained by the decay model fit).
I tried to add this by putting the following line after
https://github.com/ME-ICA/tedana/blob/c4c1621b7d0e509cd3e705771b33b2b548197468/tedana/decay.py#L539
r_squared[use_vox, :] = 1 - (np.sum((data_echo - predicted_data) ** 2, axis=1) / np.sum((data_echo - np.mean(data_echo, axis=0)) ** 2, axis=1))
Looking at the results, the R^2 values seem to be >>0.9. Since the numerator is similar to RMSE, I think this is showing that our RMSE maps might be mostly a signal dropout map and maybe we should also be saving the R^2 maps (with the variance explained separately by S0 and T2* also being saved).
I think parameter estimate variance will be more useful than the RMSE spatial map, but I'm not sure they'll be a good substitute for R^2. I've tried to calculate R^2 multiple times over the years and never gotten anything interpretable.