Icasso

[bqrosen]

I added icasso as a non-default optional method for estimating ICs in the IcaFix pipeline. After icasso is run, melodic is used for mixture modeling only and then FIX proceeds unchanged. It is currently only implemented for interpreted matlab mode. If hcp_fix_multi_run is called with --ica-method=ICASSO and {FSL_FIX_MATLAB_MODE} != 1, it will throw an error before starting the analysis.

[bqrosen]

I think I addressed most of these in a new push which also includes multilevel icasso. The matlab mode issue is thornier. The matlab script optionally includes savefig, so if that option is invoked (as it is in the hcp_fix_multi_run) then octave will fail. I’m not sure if compiled matlab can save .fig's but I think it can.

…

-Burke

On Jun 6, 2024, at 6:28 PM, Tim Coalson ***@***.***> wrote: @coalsont commented on this pull request. In ICAFIX/hcp_fix_multi_run <#294 (comment)>: > @@ -352,6 +360,10 @@ then log_Err_Abort "--fallback-threshold is not allowed when --vol-wisharts=1 is used" fi +if [[ ${FSL_FIX_MATLAB_MODE} != 1 ]] && [[ ${ICAmethod} == ICASSO ]];then Also, this would have prevented interpreted octave from trying it, so not a good test condition, even if we wanted to not recompile for some reason. — Reply to this email directly, view it on GitHub <#294 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ABPJSUYKM2FUORAHVWWMP6TZGDWBDAVCNFSM6AAAAABI5D6ULCVHI2DSMVQWIX3LMV43YUDVNRWFEZLROVSXG5CSMV3GSZLXHMZDCMBTGM2TSOBQGA>. You are receiving this because you authored the thread.

[glasserm]

niftiread is not a function I am familiar with. I typically use FSL's read_avw and save_avw to avoid adding dependencies.

[coalsont]

Technically, read_avw comes from fsl, so it was a dependency we added. niftiread is apparently part of matlab's image toolbox, starting in 2017b. Octave-forge doesn't seem to have integrated it yet, though.

[glasserm]

I didn't know that. Perhaps could check if the function is present and fall back to read_avw/save_avw if not?

[glasserm]

Is melodic_Tmodes really the same as melodic_mix?

[glasserm]

Need to do the same for the timeseries, right? PostFix may also already take care of both of these.

[glasserm]

I am assuming the contents of --mix doesn't matter, but this was my implementation to do mixture modelling on CIFTI data (perhaps that now works natively though?):

/media/myelin/brainmappers/HardDrives/2TBB/Connectome_Project/Pipelines/MixtureModelingRSN.sh

I also used a different --mmthresh=0.

[glasserm]

Please see my notes above. Also, I am not a bit fan of aggressive line continuations. They make code harder for me to read. Screens are not tiny any more...

[bqrosen]

While niftiread, niftiwrite, and niftiinfo are standard matlab toolbox functions, The script calls melodic so it’s always going to be FSL dependent. It looks like there was an attempt to add nifitread, niftiwrite, and niftiinfo to octave forge in 2022, https://savannah.gnu.org/patch/?9853, but it wan’t finished for some reason. I can swap to read_avw for better octave compatibility. That said, nifitread is at least twice as fast as read_avw for nii.gz because read_avw has to resave and reload the image to a temporary file to unzip it. nifitread can also read cifti cdata. However, niftiinfo cannot read the cifti header (it returns an intent code error). niftiwrite can write nifti1 or nifti2, compressed or uncompressed, combined or with separate .hdr.

…

-Burke

On Jun 10, 2024, at 3:14 PM, Tim Coalson ***@***.***> wrote: @coalsont commented on this pull request. In ICAFIX/scripts/run_icasso.m <#294 (comment)>: > + maxIter = str2double(maxIter); +end + +%% parse paths +% inputs +vntsFile = sprintf('%s/%s_vnts.nii.gz',ConcatFolder,concatfmrihp); +brainMaskFile = sprintf('%s/%s_brain_mask.nii.gz',ConcatFolder,concatfmri); +if ~exist(vntsFile,'file');error('%s doesn''t exist!',vntsFile);end +if ~exist(brainMaskFile,'file');error('%s doesn''t exist!',brainMaskFile);end + +% outputs +outDir = sprintf('%s/%s.ica/filtered_func_data.ica',ConcatFolder,concatfmrihp); +[~,~] = unix(['mkdir -p ' outDir]); + +%% load data, reshape to 2d, and apply brain mask +vnts = double(niftiread(vntsFile)); Technically, read_avw comes from fsl, so it was a dependency we added. niftiread is apparently part of matlab's image toolbox, starting in 2017b. Octave-forge doesn't seem to have integrated it yet, though. — Reply to this email directly, view it on GitHub <#294 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ABPJSU2VP6HB6YGZVY7DLRTZGYCKZAVCNFSM6AAAAABI5D6ULCVHI2DSMVQWIX3LMV43YUDVNRWFEZLROVSXG5CSMV3GSZLXHMZDCMBYGY3DSNJSGA>. You are receiving this because you authored the thread.

[bqrosen]

Sounds good. I’ll update the matlab code in the multiecho fitting to do that as well.

…

-Burke

On Jun 10, 2024, at 6:12 PM, glasserm ***@***.***> wrote: @glasserm commented on this pull request. In ICAFIX/scripts/run_icasso.m <#294 (comment)>: > + maxIter = str2double(maxIter); +end + +%% parse paths +% inputs +vntsFile = sprintf('%s/%s_vnts.nii.gz',ConcatFolder,concatfmrihp); +brainMaskFile = sprintf('%s/%s_brain_mask.nii.gz',ConcatFolder,concatfmri); +if ~exist(vntsFile,'file');error('%s doesn''t exist!',vntsFile);end +if ~exist(brainMaskFile,'file');error('%s doesn''t exist!',brainMaskFile);end + +% outputs +outDir = sprintf('%s/%s.ica/filtered_func_data.ica',ConcatFolder,concatfmrihp); +[~,~] = unix(['mkdir -p ' outDir]); + +%% load data, reshape to 2d, and apply brain mask +vnts = double(niftiread(vntsFile)); I didn't know that. Perhaps could check if the function is present and fall back to read_avw/save_avw if not? — Reply to this email directly, view it on GitHub <#294 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ABPJSU6ACHHWH74JXUM5EBDZGYXGBAVCNFSM6AAAAABI5D6ULCVHI2DSMVQWIX3LMV43YUDVNRWFEZLROVSXG5CSMV3GSZLXHMZDCMBYHEYDSMJVG4>. You are receiving this because you authored the thread.

[bqrosen]

I double checked and they are. I’ll make that more explicit in the code.

…

-Burke

On Jun 9, 2024, at 8:43 PM, glasserm ***@***.***> wrote: @glasserm commented on this pull request. In ICAFIX/scripts/run_icasso.m <#294 (comment)>: > +dlmwrite([outDir '/melodic_ICstats'],ICstats, 'delimiter', '\t'); + +% copy brainmask into melodic dir +copyfile(brainMaskFile,[outDir '/mask.nii.gz']) + +% save pcaD and pcaE? Looks like FIX doesn't need them + +%% calculate FTmix spectra +ts = struct(); +ts.Nsubjects = 1; +ts.ts = A_final; +[ts.NtimepointsPerSubject,ts.Nnodes] = size(ts.ts); +ts_spectra = nets_spectra_sp(ts); + +% save A_final and spectra as tab-delimited text +dlmwrite([outDir '/melodic_Tmodes'], ts.ts, 'delimiter', '\t'); Is melodic_Tmodes really the same as melodic_mix? — Reply to this email directly, view it on GitHub <#294 (review)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ABPJSU6YTLEKOKXYUNBQTTTZGUADVAVCNFSM6AAAAABI5D6ULCVHI2DSMVQWIX3LMV43YUDVNRWFEZLROVSXG5CSMV3GSZLXHMZDCMBWGQ4DKMZVGM>. You are receiving this because you authored the thread.

[bqrosen]

I checked (with diff) and changing --mix only alters the report .html’s and the log.txt so you’re right, it doesn’t matter. On --mmthresh, I think you’re right for generic mixture modeling that reporting the unthresholding z’s is probably best, according to the melodic documentation <https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/MELODIC>, 0.5 is the value to choose to emulate normal, non-mixture-modeling-only melodic usage. That is also consistent with the mailing list post <https://www.jiscmail.ac.uk/cgi-bin/webadmin?A2=fsl;6e85d498.1607> I commented. -Burke

…

On Jun 9, 2024, at 8:53 PM, glasserm ***@***.***> wrote: @glasserm commented on this pull request. In ICAFIX/scripts/run_icasso.m <#294 (comment)>: > +ts_spectra = nets_spectra_sp(ts); + +% save A_final and spectra as tab-delimited text +dlmwrite([outDir '/melodic_Tmodes'], ts.ts, 'delimiter', '\t'); +dlmwrite([outDir '/melodic_FTmix'], ts_spectra, 'delimiter', '\t'); + +fid = fopen([outDir '/components.txt'],'w'); +fprintf(fid,'%i: Signal\n',1:size(S_final,2));fclose(fid); +[~,~] = unix(['wb_command -cifti-create-scalar-series ' ... + sprintf('%s/melodic_FTmix %s/melodic_FTmix.sdseries.nii -transpose -name-file %s/components.txt -series HERTZ 0 %f',... + outDir,outDir,outDir,tr)]);% note: saving this file is the only thing tr is used for + +%% Mixture modeling +% follows recipe from https://www.jiscmail.ac.uk/cgi-bin/webadmin?A2=fsl;6e85d498.1607 +fprintf('performing melodic mixture modeling ...\n') +[~,~] = unix(sprintf('melodic -i %s/melodic_IC --ICs=%s/melodic_IC --mix=%s/melodic_mix -o %s --Oall --report -v --mmthresh=0.5',... I am assuming the contents of --mix doesn't matter, but this was my implementation to do mixture modelling on CIFTI data (perhaps that now works natively though?): /media/myelin/brainmappers/HardDrives/2TBB/Connectome_Project/Pipelines/MixtureModelingRSN.sh I also used a different --mmthresh=0. — Reply to this email directly, view it on GitHub <#294 (review)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ABPJSU6EHXIS2PD5YZQGXP3ZGUBHZAVCNFSM6AAAAABI5D6ULCVHI2DSMVQWIX3LMV43YUDVNRWFEZLROVSXG5CSMV3GSZLXHMZDCMBWGQ4TAOBUGE>. You are receiving this because you authored the thread.

[glasserm]

Does that threshold the data?

[coalsont]

This looks like it should be a bash script, I don't see any matrix math being done internally. Perhaps the only reason it is a matlab script is because the place you currently want to call it from is matlab code? Let's not use proprietary languages when there is a more obvious choice for the task.

[bqrosen]

You are right, though I would argue that it also makes it more easy to debug 🐛. It is designed to be called from matlab scripts because in most use cases (It is supposed to be a more general utility, not just a helper function for this run_icasso) its inputs will be created by a matlab script. Also, in the future an option could be added to feed it a matrix directly rather than a file path, which bash couldn't do. That said, a bash version would be good, too.

[coalsont]

If by "feed it a matrix" you mean "write that matrix to a file, run melodic, etc on it, and maybe read it back in later", I would still argue that the command sequence should be in a bash script, for more legibility, and potentially to be able to run the command sequence from a different script without having matlab/octave in the way. If you wanted it to look more matlabby, you could write a secondary matlab wrapper that writes a matrix to a file and runs the bash script. Unless a process actually uses matlab to do the math, I don't want its implementation to depend on matlab. The less matlab we use in the pipelines, the easier they are to maintain and use (compiled matlab, MCR version, etc).

[coalsont]

Again, I don't think we want this type of check in the committed version. We can compile it as a final step before merging.

[coalsont]

Using which is probably a better idea than trusting that wb_command will never treat "no arguments, so print the help info" as an error condition.

[coalsont]

I suppose -version would be fairly safe, and would also work on windows. But, as long as the melodic check uses which, might as well use the same approach for both.