diff --git a/csa_generate_figures/analyse_csa_all_models.py b/csa_generate_figures/analyse_csa_all_models.py index 07221bf2..a4932cc9 100644 --- a/csa_generate_figures/analyse_csa_all_models.py +++ b/csa_generate_figures/analyse_csa_all_models.py @@ -18,6 +18,7 @@ from charts_utils import create_experiment_folder import itertools from statsmodels.sandbox.stats.multicomp import multipletests +from statsmodels.stats.anova import AnovaRM FNAME_LOG = 'log_stats.txt' @@ -281,6 +282,13 @@ def compute_paired_t_test(data): logger.info(p_adjusted) +def compute_anova(df, depvar='std', subject='Subject', within=None, aggregate_func=None): + results = AnovaRM(data=df, depvar=depvar, subject=subject, within=within, aggregate_func=aggregate_func).fit() + p_value = results.anova_table["Pr > F"][0] + logger.info('ANOVA: {}'.format(results)) + logger.info('ANOVA p-value: {}'.format(p_value)) + + def main(): exp_folder = create_experiment_folder() print(exp_folder) @@ -434,6 +442,60 @@ def main(): annonate=annonate) logger.info(f'Number of subject in test set: {len(stds.index)}') + # Compute statistical test for soft_all CSA + ########################################### + logger.info('\nPaired ANOVA test between CSA across contrasts for soft_all') + csa_soft_all = dfs['soft_all'].copy().reset_index() + print(csa_soft_all.columns) + length = len(csa_soft_all.index.to_list()) + print('Lenght:', length) + #print(csa_soft_all) + #csa_soft_all['csa'] = csa_soft_all['DWI'] + #csa_soft_all['contrast'] = 'DWI' + oneCol = [] + subjects = [] + cols = ['dwi', 'mt-on', 'mt-off', 'T1w', 'T2star', 'T2w'] + for col in cols: + oneCol.append(csa_soft_all[col]) + subjects.append(csa_soft_all['Subject']) + csa_row = pd.concat(oneCol, ignore_index=True) + csa = pd.DataFrame() + csa['csa'] = csa_row + csa['Subject'] = np.ravel(subjects) + csa['contrast'] = 'test' + csa.loc[0:length, 'contrast'] = 'DWI' + csa.loc[length:2*length, 'contrast'] = 'MT-on' + csa.loc[2*length:3*length, 'contrast'] = 'GRE-T1w' + csa.loc[3*length:4*length, 'contrast'] = 'T1w' + csa.loc[4*length:5*length, 'contrast'] = 'T2*w' + csa.loc[5*length:6*length, 'contrast'] = 'T2w' + print(csa) + compute_anova(df=csa, depvar='csa', subject='Subject', within=['contrast']) + + logger.info('\nPaired ANOVA test between CSA across contrasts for hard_all') + csa_soft_all = dfs['hard_all'].copy().reset_index() + length = len(csa_soft_all.index.to_list()) + #print(csa_soft_all) + #csa_soft_all['csa'] = csa_soft_all['DWI'] + #csa_soft_all['contrast'] = 'DWI' + oneCol = [] + subjects = [] + cols = ['dwi', 'mt-on', 'mt-off', 'T1w', 'T2star', 'T2w'] + for col in cols: + oneCol.append(csa_soft_all[col]) + subjects.append(csa_soft_all['Subject']) + csa_row = pd.concat(oneCol, ignore_index=True) + csa = pd.DataFrame() + csa['csa'] = csa_row + csa['Subject'] = np.ravel(subjects) + csa.loc[0:length, 'contrast'] = 'DWI' + csa.loc[length:2*length, 'contrast'] = 'MT-on' + csa.loc[2*length:3*length, 'contrast'] = 'GRE-T1w' + csa.loc[3*length:4*length, 'contrast'] = 'T1w' + csa.loc[4*length:5*length, 'contrast'] = 'T2*w' + csa.loc[5*length:6*length, 'contrast'] = 'T2w' + print(csa) + compute_anova(df=csa, depvar='csa', subject='Subject', within=['contrast']) # Compare one model per contrast vs one for all ################################################ @@ -508,6 +570,5 @@ def main(): - if __name__ == "__main__": main()