Group comparison

This tutorial illustrates how to perform a statistical comparison between a group of scans and matching subjects within the normative model.

Loading and setting up subjects and normative model

First create a SOM instance with the bids directory containing your subjects of interest:

from scanometrics.core import ScanOMetrics_project
bids_directory = "~/Downloads/OASIS3"  # <-- example to be edited
SOM = ScanOMetrics_project(bids_directory)

To list the available models in the scanometrics default folder, you can run SOM.list_normative_models(). To load a given model you can run:

model_name = "Polynomial_dldirect_v1-0-3_DesikanKilliany_OASIS3_som-v0-1-0.pkl"  # <-- example to be edited
SOM.load_normative_model(model_name)

where <model_name> is the identifier for the model you want to load. Here, OASIS3 refers to the OASIS3 healthy control dataset, after removing scans that had more than 5% of abnormal regions in a Leave-one-out cross validation.

Normative models usually rely on age and other covariate values, some of which are categorical and need conversion to numerical values for their effect to be accounted for. The mapping is loaded along with the normative model, in a dictionary called cov2float. This makes sure that the mapping done when loading new subjects is the same that was used for training the normative model (e.g.: sex mapping to 0 for males and 1 for females). If your subject has new covariates compared to the training data, and which should be mapped to a numerical value, you should add the mapping function to the SOM.cov2float dictionary. For example, the sex covariate is mapped to a numerical value using the following key/value pair:

print(SOM.cov2float['sex'])
# out[]: {'M':0, 'F':1}

You are now ready to load your subjects of interest. The default is to load all subjects available in the BIDS directory (filtering of subjects to be tested can be done afterwards, see below):

SOM.load_subjects()

To load subjects data, we use the load_proc_metrics() method, which requires setting the normative model as reference for metric normalization. We use the ‘orig’ data in the normative model to do so:

SOM.load_proc_metrics(ref_metric_values=SOM.normativeModel.measured_metrics['orig'].copy(),
                      ref_metric_names=SOM.normativeModel.metric_names.copy(),
                      ref_covariate_values=SOM.normativeModel.covariate_values.copy(),
                      ref_covariate_names=SOM.normativeModel.covariate_names.copy())

Statistical comparison

To compare the loaded subjects against the normative model, you can use the test_group_differences() method:

group_analysis = SOM.test_group_differences(matching_covariates=['sex', 'scanner', 'sequence'])

As group comparisons is usually used for a certain condition, you can select the set of subjects using the group_covariate and group_label optional arguments. The group_covariate variable should be an array of integers, usually a column of SOM.covariate_values, with the same length as the number of loaded subjects. The group_label should be the value for the subjects you are interested into. For example, if you would like to test scans with severe AD in the OASIS3 dataset, you could run the following:

# Perform analysis on CDR=3 subjects (maps to 4 in cov2float)
group_analysis = SOM.test_group_differences(group_covariate=SOM.covariate_values[:, SOM.covariate_names.index('diagnosis_subtype'),
                                            group_label=4,
                                            matching_covariates=['sex', 'scanner', 'sequence'])

The group_analysis variable is a dictionary with several keys:

‘gp_residuals’: the residuals of each metric available in both the loaded data and the normative model, for each loaded subject.
‘gp-vs-gp_dfs’: the degrees of freedom for each tested metric.
‘gp-vs-gp_ts’: the value of the statistical variable for each metric.
‘gp-vs-gp_pvals’: the p-value of each test.
‘gp-vs-gp_logps’: the signed log10 of the p-value.
‘gp-vs-gp_cohen-d’: the effect size of each tested metric.

Full script

from scanometrics.core import ScanOMetrics_project
# Editable variables: bids directory, model_name, and subject_id
bids_directory = "~/Downloads/OASIS3"  # <-- to be edited
model_name = "Polynomial_dldirect_v1-0-3_DesikanKilliany_OASIS3_som-v0-1-0.pkl"  # <-- to be edited. Use SOM.list_normative_models() to list available models
group_variable = 'diagnosis_subtype'  # <-- optional, to be edited. Name of the column in SOM.covariate_values to be used for selecting subjects to analyse
group_label = 4  # <-- optional, to be edited. Selects subjects with CDR-3 as an example (corresponds to 4 in the covariate_values variable)
# Evaluation of a subject against a normative model
SOM = ScanOMetrics_project(bids_directory)
SOM.load_normative_model(model_name)
SOM.load_subjects()
SOM.load_proc_metrics(ref_metric_values=SOM.normativeModel.measured_metrics['orig'],
                      ref_metric_names=SOM.normativeModel.metric_names,
                      ref_covariate_values=SOM.normativeModel.covariate_values,
                      ref_covariate_names=SOM.normativeModel.covariate_names)
group_analysis = SOM.test_group_differences(group_covariate=group_variable,
                                            group_label=group_label,
                                            matching_covariates=['sex', 'scanner', 'sequence'])