# -*- coding: utf-8 -*- """ Created on Mon Nov 4 17:54:46 2024 @author: arefk """ import os import pandas as pd import seaborn as sns import matplotlib.pyplot as plt import statsmodels.formula.api as smf from statsmodels.stats.multicomp import MultiComparison from itertools import combinations import statsmodels import scipy from scipy.stats import dunnett # Set up warning suppression for convergence import warnings warnings.simplefilter(action='ignore') # Define user settings for the analysis selected_dialation_amount = 3 # Set dilation amount selected_qtype = "fa" # Set Qtype value selected_mask_names = ["L_cst_784_AMBA", "CRuT_MOp-RN_ipsilesional_CCcut", "CRuT_MOp-RN_contralesional_mirrored_CCcut", "CReT_MOp-TRN_ipsilesional_CCcut", "CReT_MOp-TRN_contralesional_CCcut", "CST_MOp-int-py_ipsilesional_selfdrawnROA+CCcut", "CST_MOp-int-py_contralesional_selfdrawnROA+CCcut", "TC_DORsm-SSp_ll+SSp_ul_ipsilesional_END+higherstepsize_", "TC_DORsm-SSp_ll+SSp_ul_contralesional_END+higherstepsize", "CC_MOp-MOp_cut", "OT_och-lgn_lgncut" ] # Set mask names to be used # Define simple anatomical names for the masks mask_name_mapping = { "CRuT_MOp-RN_ipsilesional_CCcut": "Rubropsinal (Ipsilesional)", "CRuT_MOp-RN_contralesional_mirrored_CCcut": "Rubropsinal (Contralesional)", "CReT_MOp-TRN_ipsilesional_CCcut": "Reticulospinal (Ipsilesional)", "CReT_MOp-TRN_contralesional_CCcut": "Reticulospinal (Contralesional)", "CST_MOp-int-py_ipsilesional_selfdrawnROA+CCcut": "Corticospinal (Ipsilesional)", "CST_MOp-int-py_contralesional_selfdrawnROA+CCcut": "Corticospinal (Contralesional)", "TC_DORsm-SSp_ll+SSp_ul_ipsilesional_END+higherstepsize_": "Thalamocortical (Ipsilesional)", "TC_DORsm-SSp_ll+SSp_ul_contralesional_END+higherstepsize": "Thalamocortical (Contralesional)", "CC_MOp-MOp_cut": "Corpus Callosum", "OT_och-lgn_lgncut": "Optic" } # Define the PlotBoxplot function def PlotBoxplot(filtered_df, qq, mm, dd, output_folder, sig): """ Plot boxplots for filtered dataframe and save as images. Parameters: filtered_df (DataFrame): Filtered dataframe based on "Qtype" and "mask_name". qq (str): Qtype value. mm (str): Mask_name value. output_folder (str): Path to the output folder to save images. """ # Set the color palette to grey sns.set_palette("rocket") # Set the figure size to 18 cm in width plt.figure(figsize=(18/2.54, 4)) # Convert 18 cm to inches # Create a boxplot sns.boxplot(data=filtered_df, hue="Group", y="Value", x="merged_timepoint", fliersize=3, flierprops={"marker": "o"}, palette="rocket") # Set title and labels mask_label = mask_name_mapping.get(mm, mm) plt.title(f'{qq}_{mask_label}_{dd}') plt.xlabel('Timepoint') plt.ylabel(qq + "(n.a)") # Show legend without title and frame plt.legend(title=None, frameon=False, bbox_to_anchor=(1.05, 1), loc='best') plt.tight_layout() # Save the plot as an image if sig: output_path = os.path.join(output_folder, f'{qq}_{mm}_{dd}_mixedlm_significant.png') else: output_path = os.path.join(fa_over_time_plots_dir, f'{mm}.png') plt.savefig(output_path, dpi=100) # Close the plot to avoid displaying it plt.show() plt.close() # Get the directory where the code file is located code_dir = os.path.dirname(os.path.abspath(__file__)) # Get the parent directory of the code directory parent_dir = os.path.dirname(code_dir) # Define the path for the input CSV file input_file_path = os.path.join(parent_dir, 'output', "Quantitative_outputs","onlyCST", 'Quantitative_results_from_dwi_processing_withOutWM_and_partialCST.csv') # Define the path for the output folder to save plots plots_output_dir = os.path.join(parent_dir, 'output', 'Figures', 'pythonFigs') fa_over_time_plots_dir = os.path.join(parent_dir, 'output', 'Figures', 'fa_over_time_plots') os.makedirs(fa_over_time_plots_dir, exist_ok=True) os.makedirs(plots_output_dir, exist_ok=True) # Load the CSV file into a dataframe df = pd.read_csv(input_file_path) #df = df[df["cluster ii"]=="High ii"] df = df.drop_duplicates() df = df[df["merged_timepoint"].isin([0,3,7,14,21, 28])] # Create lists to store results and terminal output results = [] # Iterate over the user-defined "dialation_amount", "Qtype", and "mask_name" for dd in [selected_dialation_amount]: for qq in [selected_qtype]: for mm in selected_mask_names: # Filter the dataframe based on current "Qtype", "mask_name", and "dialation_amount" filtered_df = df[(df["Qtype"] == qq) & (df["mask_name"] == mm) & (df["dialation_amount"] == dd)] # Check if the filtered dataframe is not empty if not filtered_df.empty: # Plot the boxplot for the filtered dataframe # Assuming `sig` is a boolean that indicates statistical significance (set to False for now) sig = False PlotBoxplot(filtered_df, qq, mm, dd, plots_output_dir, sig) # Print the list of results (currently empty as calculations are not included) print(results) # Replace this with actual results dataframe if available