Loading an image as a statistic_image and then trying to write it back to file in MNI space fails with the following error:

A simple work-around is to load it as an fmri_data object, and then write to disk; this works.

This error is likely limited to using the 'mni' flag w/ the statistic_image object

Hi,

A quick question about design matrix: when running the GLM, is it better to concatenate the data (various runs), the movement regressors, the wm and csf avg. signal in the matrix?

In some but not all cases, resample_space causes an image to be shifted left/right / in the axial plane by a few mm. I have not figured out the boundary conditions of when this happens and when it does not happen. Here is an example of when it happens:

The top is the original image. The bottom is resampled to a new space -- notice how the image has been shifted right. In my testing, this occurrence of this bug may depend more on the target space being sampled to, rather than the resampling algorithm, but more testing is needed to better understand this.

Here is the code to recreate the bug:

limbicmask = fmri_data(fullfile(repodir, 'data', 'Etienne2016modules.nii'))
limbic_rs = resample_space(limbicmask, gbt, 'nearest');
orthviews_multiple_objs({limbicmask, limbic_rs})

Contact me directly for the brain imaging files used to above to reproduce the bug; I don't have permission to post them publicly. Using other files, the bug may or may not appear.

At line 51, significant voxels are identified as

z(z < u) = 0

In the case of two-tailed tests, this eliminates negative activations. I can fix this issue by changing it to

z(abs(z) < u) = 0

but am not sure that won't break other interdependencies.

Extracting parcel local signature patterns

NPS
Arrays have incompatible sizes for this operation.

Error in image_vector/apply_parcellation (line 260)
wh = wh_in_dat & ~parcels.removed_voxels;

Error in fmri_data/extract_measures_batch (line 361)
local_pattern = apply_parcellation(data_obj_atlas_space, parcel_obj, 'pattern_expression', sig_obj);

Hi,
I have an fMRI dataset that I preprocessed using SPM12. Now, I would like to perform multilevel glm with the 'glmfit_multilevel.m' script. However, when I use the example data (as defined in lines 69-81 of the glmfit_multilevel.m), I receive the following errors:

Undefined operator '*' for input arguments of type 'cell'.

Error in glmfit_general>get_weights_based_on_varcomponents (line 506)
SSW = dfwithin' * s2within; % 1 x nvars, within-subjects sum sq. errors. dfwithin must be row vector

Error in glmfit_general (line 147)
[W, btwn_var_est] = get_weights_based_on_varcomponents(s2within, s2between_ols, dfwithin, verbose);

Error in glmfit_multilevel (line 232)
stats = glmfit_general( ...

Thanks very much in advance.

Ezgi

atl = load_atlas('canlab2023_fine_fmriprep20')
r = table(t,  'atlas_obj', atl);

Produces the following output after the results table:

Regions labeled by reference atlas CANLab2023_fine_MNI152NLin2009cAsym_2mm
                                                                          
Volume: Volume of contiguous region in cubic mm.
MaxZ: Signed max over T
Atlas_regions_covered: Number of reference atlas regions covered at least 25% by the region. This relates to whether the region covers 
multiple reference atlas regions                                                                                                       
Region: Best reference atlas label, defined as reference region with highest number of in-region voxels. Regions covering >25% of >5 
regions labeled as "Multiple regions"                                                                                                
Perc_covered_by_label: Percentage of the region covered by the label.
Ref_region_perc: Percentage of the label region within the target region.
modal_atlas_index: Index number of label region in reference atlas
all_regions_covered: All regions covered >5% in descending order of importance
 
For example, if a region is labeled 'TE1a' and Perc_covered_by_label = 8, Ref_region_perc = 38, and Atlas_regions_covered = 17, this means 
that 8% of the region's voxels are labeled TE1a, which is the highest percentage among reference label regions. 38% of the region TE1a is  
covered by the region. However, the region covers at least 25% of 17 distinct labeled reference regions.                                   
                                                                                                                                           
References for atlases:
                       
 
&
(
)
,
-
.
/
0
1
2
3
4
5
6
7
8
9
:
;
?
A
B
C
D
E
F
G
H
I
J
K
L
M
N
O
P
R
S
T
U
V
W
X
Y
Z
[
]
a
b
c
d
e
f
g
h
i
j
k
l
m
n
o
p
r
s
t
u
v
w
x
y
z
é
í
ö
–
                                                                                                     
Note: Region object r(i).title contains full list of reference atlas regions covered by each cluster.

Writing them to disk with the 'mni' flag creates rings of zero and near-zero values around positive activation clusters. Here is the same image written to disk with and without the 'mni' flag, and then re-loaded:

To reproduce: I'm attaching a sample file that can be loaded in as a statistic_image, written to disk with and without the 'mni' flag, and then reloaded
tspecific.nii.gz

I've found a weird memory leak in the roc_plot() function when running in single-interval mode and the accuracy is perfect. Don't have time to figure it out now, but will try and come back later.

binary_outcome = logical([ones(sum(dat.Y==1),1); zeros(sum(dat.Y==3),1)]);
pattern_response = [pexp(dat.Y==3); pexp(dat.Y==1)];

% Single Interval Accuracy
ROC_si = roc_plot(pattern_response, binary_outcome, 'color', 'r');

Hi there,

I’m trying to verify that the results I get from glmfit_multilevel and that I get from lme4 and nlme in R are the same, but I am finding that they are not. My model has one continuous predictor and one categorical predictor, both on level 1, and I can include a continuous interaction as well. I’m new to using glmfit_multilevel, but have been using R for years, so am wondering how glmfit_multilevel computes betas and SE that is different from the R packages. I’ve attached a reproducible example. Any clarification on this would help.

Thanks in advance,

Liz

Glmfit_multilevelvsR.zip

To reproduce error:

a = [
19.0000 29.6667 21.3333 23.4667
14.6667 27.9333 53.0667 60.6667
58.9333 61.4667 70.8000 75.6667
70.8667 85.5333 75.6667 89.6667
25.0667 33.2000 51.0000 74.2000
42.8667 65.9333 76.6000 93.4000
19.2000 31.4667 55.4667 55.6000];

barplot_columns(a, 'noviolin');

OUTPUT:

Index exceeds matrix dimensions.

Error in barplot_columns (line 515)
mycolcolor = mycolor(i, :) ./ 2;

using the 'cmaprange' option in addblobs.m gives a warning;

'Warning: Unknown input string option:cmaprange '

and blobs of some maps are completely black when using the cmaprange option now. Works for some maps normally, although still issues a warning.

I have experienced erratic behavior with the legend function when used for statistic_image objects, in terms of layout, anchors, and colormaps, whether standalone or as called by the montage function.

Has anyone experienced the same problems and how could we fix this?

If two atlas being merged contain an overlapping set of regions, the merged atlas contains double copies of this region.

Toy example that reproduces the bug:

atlas_obj = load_atlas('canlab2018_2mm');
amy = select_atlas_subset(atlas_obj, {'Amy'})
double_amy = merge_atlases(amy, amy)

% check result
double_amy.labels

ans =

  1×8 cell array

    {'Amygdala_CM_'}    {'Amygdala_SF_'}    {'Amygdala_AStr_'}    {'Amygdala_LB_'}    {'Amygdala_CM_'}    {'Amygdala_SF_'}    {'Amygdala_AStr_'}    {'Amygdala_LB_'}

The 'always_replace' and the 'noreplace' options to merge_atlases do not address this.

I hesitate to try and fix because there already seems to be a complex logic regarding the 'always_replace' and the 'noreplace' options, and I am unsure how this would intersect with the issue reported here.

Dear Tor and co-developers of Mediation Toolbox,

I am contacting you to ask a few questions rather than raising an issue with the toolbox :) I am using the mediation toolbox for analysing a longitudinal (2 sessions - 3 groups) multimodal (functional, structural, DTI) brain imaging+behavioural data. The data has been preprocessed, therefore, we are looking into using the toolbox for not the raw data (as per examples provided with the toolbox) but the preprocessed one. We have multiple potential mediators (large size data - 67 patients x 87 values or features) and I am concatenating the matrices into one large variable (e.g., M = [Vol1 Vol2 fc1 fc2 sc1 sc2 ...], 67 x N). However, I have a few questions:

1. The toolbox seems to be able to deal with categorical/multinomial (either binomial or more) variables as the independent variable (X - intervention/patient groups) and also as covariates (sex, centre, etc.). However, we seem to have issues/errors with this... of course, the two variables for M and Y have been standardised.
2. We are very interested in finding out which columns of the matrix M (mediator as explained above) are chosen as potential mediators and how good of mediators they are individually. However, we just have two variables, being paths and stats, as outputs. Is there anyway of finding out how the features/columns of matrix M perform as mediators and which ones are the best ones?
3. Unfortunately, we do not seem to be able to plot the outputs either. This will help with better understanding the underlying elements of the mediation process.

At the moment, the piece of code that we are using is as follows:

[paths, stats] = mediation(X, Y, M, 'boot', 'plot', 'verbose', 'covs', covariates, 'covnames', covnames, 'bootsamples', 10000);

Any suggestions on how we can improve our analysis or find out more about the details of the mediation would be very much appreciated.

Warmest wishes and please stay safe,
Amir Dehsarvi
Research Fellow | Aberdeen Biomedical Imaging Centre | School of Medicine, Medical Sciences, and Nutrition | University of Aberdeen

With the 'full' montage type, removeblobs() will remove blobs from all the other montages, but not from the surface displays (the top row)

Dear Canlab group, thank you for this great toolbox! I am using your scn_session_spike_id function and it works very nicely.
However, using the default MADs value of 10, quite many of my volumes are identified as containing spikes (up to 10 % from one session). This made me wonder what your rationale was for setting this default threshold? I first thought that I could simply find the correct threshold by determining the critical chi-square value for a desired probability level and the degrees of freedom. Then I saw that in your function the MADs threshold is not applied directly to the calculated mahalanobis distances in d2 (I guess given the degrees of freedom it would also have to be much higher than 10 in this case anyways) but rather as a cutoff for deviation of these d2 values from a moving average. Can you give a recommendation what criteria to use if one wants to change the default MADs value?

Many thanks in advance!

Teresa

I am configuring a new computer and can't seem to find the brainmask.nii file, which is the default mask used by fmri_data(). Should we put this somewhere in the repository (e.g., canlab_canonical_brains or External)? or put a note on where to get it? I have spm12 and spm8 and it doesn't seem to be in the canonical folder anymore.

I'm going to put it in canlab_canonical_brains for now, but feel free to move it.

Hi,

Statistic_images crash when attempting to resample_space():

This is in fact a semi-known issue, at least to Wani and me (but has not been officially logged). For example, image_vector.write() will resample the image space when the MNI flag is passed in. And, if the input to write() is a statistic_image, we first cast it to fmri_data before resampling the space -- because resampling the space of a statistic_image crashes.

Also: I think it would be super helpful to have a brief write up of how we are storing image space parameters (i.e., a description of volInfo, the image "mask", etc.) -- and why statistic_image does this differently, apparently.

Thanks
Yoni

Hi

I'm wondering if there are any legend information regarding the 489 regions in the combined 2018 atlas. So I want to use this atlas in FSL however the index are numbers 1-489. Are there any corresponding labels somewhere? Thank you.

Kevin

Hi CanlabCore maintainers,

I have a question about the calculation of list lengths for fMRI stimulus sequences in OptimizeDesign11/optimizeGA. I can see in line 233 that when restevery and restlength are turned on in the controller script (and dorests is positive), numRestStim is calculated as the total number of ISI-length rest intervals to add onto the duration of the task (including stimulus intervals and blank ITI intervals).

My question is, what function is that - 1 performing in that formula? I was using optimizeGA to optimize some run orders recently, and I noticed that the numbers of trials and rest intervals that get reported out by the optimizeGA output weren't lining up exactly with what I expected based on what I input for freqConditions, scanLength, restevery, and restlength. (I was calculating scanLength as a function of condition counts and rest intervals, so it should have lined up exactly!) However, in our local version of the CanlabCore repository, I tried removing the - 1 from the formula calculating numStim with rests, and then the numbers coming out of optimizeGA lined up with the numbers going in. (The break intervals that show up in the resulting optimized trial orders were the right length as well after I removed the - 1.)

I wanted to write here to make sure I wasn't breaking something else in the formula? Thanks for your guidance!

On such regions such as NPS 'rLOC' and SIIPS 'L-TP', atlas2region () doesn't work for certain regions.

Example:

Step 1: Download NPSneg and SIIPSneg from this link: https://www.notion.so/Generate-and-write-out-NPS-and-SIIPS-ROI-Subregion-Objects-8eb8e2f3bee34dd2bf8bc0477ef6e3ec

Step 2:
load('npsneg.mat');
load('siipsneg.mat');
at=npsneg;
atlas2region(at.select_atlas_subset({'rLOC'}))
at=siipsneg;
atlas2region(at.select_atlas_subset({'L-TP'}))

When applying pain_cog_emo masks, some signatures always yield 0s. I traced it back to an error in a few lines in the code.

Lines 195 ... 197 in apply_mask.m:
to_remove = ~inboth(dat.volInfo.wh_inmask);
to_remove_mask = ~inboth(mask.volInfo.wh_inmask);

wh_inmask is always a single column
inboth is multi-column when applying multiple masks and only its first column is included in the current code.

I tried replacing line 171 from
inmaskdat = logical(mask.dat);
to
inmaskdat = get_nonempty_voxels(mask);
but this leads to improper norms for the data the mask is applied on.

Maybe the best way would be to check if multiple masks and run apply_mask again for each individually.

Using 'rbf' flag in 'cv_svm' for predict, the predicted values (stats.yfit) is all 1's. That is, it predicts all observations as belonging to ONE class.

Without the radial basis function, some observations are predicted to be class A, others as class B, as expected.

The file for reproducing the error is in the repo, see CanlabPrivate/FILES_FOR_REPRODUCING_BUGS/svm_rbf_bug

Hi,

Nikhil and Tor -- thank you for the new display style functionality!!
When I make a "full" montage, and then add blobs to it from a statistic image, they only show on the axial montage, not on the sagital/surface montages at the top row. This is true whether I try to add them using addblobs() or with canlab_results_fmri_display()

Hi
I have a question regarding mediation toolbox. I want to see whether time series related to region "M" in the brain is mediating the relation between time-series related to region "X" and region "Y"! Actually, I am not sure which preprocesing steps are necessary! In my case, Time-Series are scaled to percent signal change. In addition I have used detrending, mean centering, and despiking before running mediation function. The problem is that the values for X->M and M->Y and X->Y are not in line with correlation values between related signals! I understand that they are coefficients for multiple regression analysis but again they should be similar to correlation results! Therefore I just want to make sure which preprocessing steps are essential.
Thank you for your time.
Best Regards,
Neda

Hi ,

I used mediation_brain_results for Fdr analysis, and I noticed that in the output I get a warning as:
'p-values in cl.Z will not give valid spm_max subclusters'.

t-values, and z values are also not reported for the peaks in the outputs of FDR, and only neg(logp) is reported as maxstats.

Does this warning implies that the neg(logp) reported as maxstats is not accurate?
Is the table reported with the coordinates of voxels that survive FDR, valid?

I also run the analysis with mediation_brain_results and cluster thresholding, and in the reports of clusters and subclusters, I have the t value, z value, and neglogp for all the peaks. I was just wondering if the t values reported for subclusters are accurate?

some of the peaks in FDR matches with the peaks in subclusters and they both have the same neg(logP), so I assume that in the table report
of subclusters, the neg(logP) is not accurate, as well.

Thanks for clarifying these few issues,

Mina

If DSGN.allowmissingfunc = true, and you have missing functional files for runs that are not the last runs of your subject, noise_regs.mat files (containing the spm-style R matrix) are not assigned to the correct runs.
When a folder for a run is empty, this is handled correctly for onsets and durations files, but not for the noise_reg files. For the latter, the code looks for a noise_regs.mat file in the folder for the empty run, does not find one, causing the mutliple regressors field in the spm batch file to be empty, and then picking the wrong noise_regs.mat file for the following run (the file for run n-1 is assigned to run n).
The problem seems to be located in canlab_glm_subject_levels_run1subject.m, which is a child process of canlab_glm_subject_levels.m - I may need some help to fix it but should definitely be feasible @bogpetre @torwager

Hi,

Can anyone explain to me the error check in region.m, line 265, pasted below? The code below the error check seems to work fine, even when the error check fails. I also don't understand conceptually the error check. Maybe there is an error in the error check?

thanks, Y

        % If statistic_image, we need to consider thresholding
        % so apply .sig field
        if isa(mask, 'statistic_image')                
            if size(mask.dat, 1) == mask.volInfo.nvox
                error('statistic_image objects should not have data vector (.dat) the length of full image space.')
            end

            maskData = maskData .* mask.sig(:, 1);

        end

Upon running the canlab_toolbox_setup I get this warning

Warning: Function vecnorm has the same name as a MATLAB builtin. We suggest you rename the function to avoid a
potential name conflict. 
> In path (line 109)
  In addpath (line 96)
  In canlab_clone_github_repository (line 164)
  In canlab_toolbox_setup (line 179) 

Might be worth renaming it to avoid name-space collision.

When balanced accuracy is requested in the function roc_plot, the output p values (ROC.accuracy_p) do not account for difference in class frequencies (the standard binomial test against .5 is computed).

This leads to cases where balanced accuracy is at .5, but p-values can be highly significant. An alternate method for computing p should probably be implemented when balanced accuracy is requested.

running addbrain('surface left) or addbrain('surface right') returns:
I need the file: L.pial_MSMAll_2_d41_WRN_DeDrift.32k.mat

From the documentation, this is the template "from Glasser et al. 2016 Nature based on the Human Connectome Project." Is it intended to be included in CanlabCore?

Dear Canlab team,

When I try to run the canlab_glm_subject_levels function I get the following error: 'undefined function 'canlab_spm_fmri_model_job' for input arguments of type 'cell''.

The problem seems to be that this function is missing (from the matlab path). However, upon checking I cannot find this function (as well as the function 'canlab_spm_contrast_job_luka') in any of the canlab repositories. Is this correct?

I added a screenshot of the DSGN.mat variables in case the error may come from there.

Many thanks in advance,

Jelle

Hello Canlab,

I'm looking forward to using core tools for my data but have searched through the canlab repository and can't find any of the preproc scripts. Have the names been changed or have they not been added to the repository? I use spm12 or 8. I appreciate your time!

When the number of hits is equal to exactly half the number of trials, roc_plot returns a p-value greater than 1.

The error is in the function binotest.m line 24. To simulate, simply run:

binotest([0 0 1 1]', .5)

Its taking the cdf of getting at least 1/2 hits, which is greater than .5, and then multiplying that value by 2, yielding a pval > 1. I'm not sure what the right fix is.

Hello Tor and co-developers,

I am new to running mediation_brain and encountered an issue.
While running mediation_brain, we intend to use brain images as the mediator with X independent variable and Y dependent outcome score.
In the run at the stage of mediator search, it reads the brain images list provided as mediators and creates an array with the values of every voxel in it.
In the mediation part, inputs get set up that matches the length of all three inputs X, Y, and M. At this point, the error occurs as the unequal length of the variables X, Y and M the following:
"Error using mediation>setup_inputs (line 1223)
X, Y, and M must be equal in length".

Any suggestions on how we can resolve the issue or find out more about the details of the mediation would be much appreciated.

Thank you

Apurva

icatb_fastICA is a dependency not checked for the image_vector.ica() method. icatb_fastICA is found in the GIFT package in MatLab, but there needs to be an easier way to inform users of this dependency.

movie_tools and movie_of_slice_timeseries are no longer functional because avifile() has been replaced by writeVideo().

timeseries_extract_slice, which is called in movie_of_slice_timeseries, returns all 0s.

Labeling errors were found on CANlab_combined_atlas_object_2018_2mm.mat. In modeling_labeling_descriptions, cortex_visual_central and cortex_visual_peripherial are inverted, and need to be modified. I attached a fixed atlas that I modified.

CANlab_combined_atlas_object_2018_2mm.mat.zip

Dear Canlab team,

since I've been working with this script recently, I just wanted to point out to you that the current version of the function onsets2fmridesign doesn't work (at least in recent Matlab versions, I'm using Matlab R2019b). The reason seems to be a change in the function pad which has been changed to padarray (see also here: fordanic/image-registration#6)

If I replace 'pad' with 'padarray' in line 456, the function works just fine. I think this could be helpful for future users.

Thanks.

Best wishes
Sandra

Hello there,

I was recently trying a glmfit. I tested the example given in the help for the function (glmfit_multilevel) and I get the following error. It seems that the error is rather deep in the code and I am not good enough with statistics to determine how to fix this. Has something changed since the example was written and tested?

Undefined operator '*' for input arguments of type 'cell'.

Error in glmfit_general>get_weights_based_on_varcomponents (line 506)
SSW = dfwithin' * s2within; % 1 x nvars, within-subjects sum sq. errors. dfwithin must be row vector

Error in glmfit_general (line 147)
[W, btwn_var_est] = get_weights_based_on_varcomponents(s2within, s2between_ols, dfwithin, verbose);

Error in glmfit_multilevel (line 232)
stats = glmfit_general( ...

Many thanks,
Rob

When using the 'cv_svr' option in fmri_data predict, the predicted intercept seems to always be 0. This is true for the full-data pattern and for each CV fold's pattern. In the output, see the 2nd column in stats.other_output and in stats.other_output_cv

Resample space incorrectly transforms data if written in any order other than x, y, z dimensions. The transformation assumes this order, but it is not necessarily the case.

Hi everyone,
I would like to run 'scn_spm_design_check.m', but I am getting a lot of error like 'Undefined function or variable 'getvif''. Is it because of the version of SPM(12) or Matlab (17b) I use?

Thanks
Thomas

I have canlabcore and the Robust toolbox and kept getting an error that it could not find the CreateEXPT function.

I fixed it by searching and digging that out of another file called SCN_Core_Support found on a different github, but it seems to duplicate most of what is in CanLabCore.

Should EXPT_setup_tools be included in CanLabCore?

Dear Canlab,

I'm trying to load firstlevel zmaps into a canlab dataobj to use the multivariate prediction tools.

When using:

maps = dir('/data04/Sebastian/RewardSignature/analysis/Results/Univariate/Firstlevel/HR-won/concat/*res*.gz') data_obj = fmri_data(maps);

I'm getting the following error:

Reference to non-existent field 'volInfo'.

Error in fmri_data (line 442)
obj.mask.volInfo = obj2.volInfo;

when checking the header of my zmaps there is information on the pixel dimensions as well as the dimensions of the data.

the fmri_data() call works fine when running on prerocessed data. So the issue seems to be specific to loading statistical maps.

Is it possible at all to load these maps into a dataobj?

Thanks for your help,

Sebastian

It was nice for users outside of the laboratory to be able to use our toolboxes since there is no documentation about dependencies anywhere. I had checked the licenses for all of the toolboxes that were in the external folder and it seemed like we were allowed to redistribute as long as it was noted who developed the toolbox and where we acquired it.

I don't really care either way, but someone should create documentation for installing the CanlabCore repository and include a list of dependent toolboxes and where they can be downloaded.