<jupyter_start><jupyter_text>Permutation of group differences.<jupyter_code>import os
import pandas as pd
import numpy as np
import networkx as nx
import scona as scn
from scona.visualisations import plot_network_measures, plot_degree_dist, plot_rich_club
import matplotlib.pyplot as plt
import nilearn.plotting as plot
import functions.plotting_functions as Pfun
import functions.statistical_functions as Sfun
import seaborn as sns
sns.set_style('dark')<jupyter_output><empty_output><jupyter_text>This code creates the enviornmental variable for where the data is stored. Create a .env file in the directory with the file path to data.<jupyter_code>from decouple import config

data = config('data')
lh_volume = pd.read_csv(f'{data}/lh_volume.dat',sep='\t').drop([
                     'BrainSegVolNotVent', 'eTIV'],axis=1).rename(columns={'lh.aparc.volume':'G-Number'})

rh_volume =  pd.read_csv(f'{data}/rh_volume.dat',sep='\t').drop([
                       'BrainSegVolNotVent', 'eTIV','rh.aparc.volume'],axis=1)

group = pd.read_csv(f'{data}/cortical_measures.csv').iloc[0:,2]

volume = pd.concat([lh_volume, rh_volume, group],axis=1)

names = list(volume.columns.drop(['G-Number','age_adjusted_group']))

centroids = pd.read_csv(f'{data}/atlas.csv')

centroids = centroids[['x.mni', 'y.mni', 'z.mni']].to_numpy()

group = volume.groupby('age_adjusted_group')
aan = group.get_group('AAN').reset_index(drop=True)
hc = group.get_group('HC').reset_index(drop=True)
wr = group.get_group('WR').reset_index(drop=True)
aan.head()
aan_graphs = Sfun.create_graphs(aan.iloc[:,1:69],names, centroids)
wr_graphs = Sfun.create_graphs(wr.iloc[:,1:69],names, centroids)
hc_graphs = Sfun.create_graphs(hc.iloc[:,1:69],names, centroids)
brain_bundle_aan = scn.GraphBundle([aan_graph_threshold], ['AAN_graph_thresholded'])
brain_bundle_aan.create_random_graphs('AAN_graph_thresholded', 1000)
brain_bundle_wr = scn.GraphBundle([wr_graph_threshold], ['WR_graph_thresholded'])
brain_bundle_wr.create_random_graphs('WR_graph_thresholded', 1000)
brain_bundle_hc = scn.GraphBundle([aan_graph_threshold], ['HC_graph_thresholded'])
brain_bundle_hc.create_random_graphs('HC_graph_thresholded', 1000)
aan_global_measures = brain_bundle_aan.report_global_measures() 
aan_global_measures
hc_global_measures = brain_bundle_hc.report_global_measures() 
hc_global_measures
wr_global_measures = brain_bundle_wr.report_global_measures() 
wr_global_measures
Pfun.distro_plots(aan_global_measures)
measures = ['average_clustering', 'average_shortest_path_length', 'assortativity', 'modularity', 'efficiency'] 
measure = measures[1]
tcrit_value = aan_global_measures[measure].mean() - hc_global_measures[measure].mean() #- wr_global_measures[measure].mean()


permutations = 10000

combined = np.concatenate([aan_global_measures[measure].values, hc_global_measures[measure].values]) #wr_global_measures[measure].values, 

index_1 = len(aan_global_measures[measure].values)
#index_2 = len(wr_global_measures[measure].values) + len(wr_global_measures[measure].values)

ctvalue=[]
for perm in range(permutations):
        np.random.shuffle(combined)
        random_group_one = combined[0: index_1]
        #random_group_two = combined[index_1: index_2]
        random_group_three = combined[index_1:]
        crit_val = random_group_one.mean() - random_group_three.mean() #- random_group_two.mean() 
        ctvalue.append(crit_val)

array = np.array(ctvalue)
mean = array.mean()
std_dev = np.std(ctvalue)
lower_2std = mean - 2*std_dev
upper_2std =  mean + 2*std_dev
lower_1std = mean - std_dev
upper_1std = mean + std_dev

if tcrit_value >= mean:
        sum = [val for val in ctvalue if val >= tcrit_value]
        print('upper')
else:
        sum = [val for val in ctvalue if val <= tcrit_value]
        print('lower')
     
     #if tcrit_value > 0 and mean > 0:
     #           sum = [val for val in ctvalue if abs(val) <= abs(tcrit_value)]
     #           print('lower 1')
     #else:
     #           sum = [val for val in ctvalue if abs(val) >= abs(tcrit_value)]
     #           print('lower 2')

        
        
pval = len(sum) / permutations



print(pval)

print(tcrit_value)
print([lower_2std, upper_2std])
hist = sns.displot(data=ctvalue, height=10).set(title=f'Histplot for {measure} after {permutations} permuatations')
hist.refline(x=tcrit_value, color='purple')
hist.refline(x=lower_2std, linestyle='-', color='black')
hist.refline(x=upper_2std, linestyle='-', color='black')
hist.refline(x=lower_1std, linestyle='-', color='red')
hist.refline(x=upper_1std, linestyle='-', color='red')
plt.show()
measures = ['average_clustering', 'average_shortest_path_length', 'assortativity', 'modularity', 'efficiency'] 
measure = measures[4]

tcrit_value = aan_global_measures[measure].mean() - wr_global_measures[measure].mean() - hc_global_measures[measure].mean() 

permutations = 10000

combined = np.concatenate([aan_global_measures[measure].values, wr_global_measures[measure].values, hc_global_measures[measure].values]) 

index_1 = len(aan_global_measures[measure].values)
index_2 = len(wr_global_measures[measure].values) + len(wr_global_measures[measure].values)

ctvalue=[]

for perm in range(permutations):
        np.random.shuffle(combined)
        random_group_one = combined[0: index_1]
        random_group_two = combined[index_1: index_2]
        random_group_three = combined[index_1:]
        crit_val = random_group_one.mean() - random_group_two.mean() - random_group_three.mean()  
        ctvalue.append(crit_val)

mean_val = Sfun.mean_std(ctvalue)


if tcrit_value >= mean_val['mean']:
        sum = [val for val in ctvalue if val >= tcrit_value]
        print('upper')
else:
        sum = [val for val in ctvalue if val <= tcrit_value]
        print('lower')
     
     #if tcrit_value > 0 and mean > 0:
     #           sum = [val for val in ctvalue if abs(val) <= abs(tcrit_value)]
     #           print('lower 1')
     #else:
     #           sum = [val for val in ctvalue if abs(val) >= abs(tcrit_value)]
     #           print('lower 2')

        
        
pval = len(sum) / permutations



print(pval)

print(abs(tcrit_value))
print(abs(mean_val['mean']))
Pfun.pval_plotting(ctvalue, tcrit_value, measure, permutations)

def monte_carlo(group_one, group_two, group_three, permutations, plotting=False, measure='None'):
    
    '''
    
    Function to run a monte_carlo permuation hypothesis testing.
    
    Calculates a critical value then randomly shuffles groups calculating the crticial value 
    for each permutation.

    Pval calculated as:
    sum of all values greater than critical value / number of permutations

    Parameters
    ----------------------------------------------
    group one: numpy array of values for group one
    group one: numpy array of values for group two
    group one: numpy array of values for group three
    permutations: int, number of permutations to perform.

    Returns
    ---------------------------------------------
    results: dict object with critical value, pvalue and std
    
    '''

    critical_value = group_one.mean() - group_two.mean() - group_three.mean()
    combined = np.concatenate([group_one, group_two, group_three])

    index_1 = len(group_one)
    index_2 = len(group_two) + len(group_two)

    ctvalue = []

    for perm in range(permutations):
        np.random.shuffle(combined)
        random_group_one = combined[0: index_1]
        random_group_two = combined[index_1: index_2]
        random_group_three = combined[index_2:]
        crit_val = random_group_one.mean() - random_group_two.mean() - random_group_three.mean()
        ctvalue.append(crit_val)
    
    array = np.array(ctvalue)
    mean = array.mean()
    std_dev = np.std(ctvalue)
    lower_2std = mean - 2*std_dev
    upper_2std =  mean + 2*std_dev

    if tcrit_value >= mean:
            sum = [val for val in ctvalue if val >= tcrit_value]

    else:
            sum = [val for val in ctvalue if val <= tcrit_value]
             
    pval = len(sum) / permutations



print(pval)

print(abs(tcrit_value))
print(abs(mean))
#print([lower_2std, upper_2std])
hist = sns.displot(data=ctvalue, height=10).set(title=f'Histplot for {measure} after {permutations} permuatations')
hist.refline(x=tcrit_value, color='purple')
hist.refline(x=lower_2std, linestyle='-', color='black')
hist.refline(x=upper_2std, linestyle='-', color='black')
plt.show()

    result = {'critval':critical_value, 'pval': pval,'lower_std':lower_2std, 'upper_std':upper_2std}

    if plotting == True:
        pval_plotting(ctvalue, critical_value, measure, permutations)

    return result
measure = 'assortativity' 
plots_for_pval = monte_carlo(aan_global_measures[measure].values, wr_global_measures[measure].values, hc_global_measures[measure].values, permutations=50000, plotting=True, measure=measure)
for measure in aan_global_measures.columns: 
    eff = monte_carlo(aan_global_measures[measure].values, wr_global_measures[measure].values, hc_global_measures[measure].values, permutations=50000)
    print(eff)
def monte_carlo_2(group_one, group_two, permutations):
    
    '''
    Function to run a monte_carlo permuation hypothesis testing.
    
    Calculates a critical value then randomly shuffles groups calculating the crticial value 
    for each permutation.

    Pval calculated as:
    sum of all values greater than critical value / number of permutations

    Parameters
    ----------------------------------------------
    group one: numpy array of values for group one
    group one: numpy array of values for group two
    group one: numpy array of values for group three
    permutations: int, number of permutations to perform.

    Returns
    ---------------------------------------------
    results: dict object with critical value and pvalue
    
    '''

    critical_value = group_one.mean() - group_two.mean()
    combined = np.concatenate([group_one, group_two])

    ctvalue = []

    index_1 = len(group_one)

    for perm in range(permutations):
        np.random.shuffle(combined)
        random_group_one = combined[:index_1]
        random_group_two = combined[index_1:]
        crit_val = random_group_one.mean() - random_group_two.mean() 
        ctvalue.append(crit_val)
    
    sum = [val for val in ctvalue if abs(val) >= abs(critical_value)]
    
    #print(len(sum))
    #print('\n', ctvalue,'\n')

    pval = len(sum) / permutations

    result = {'critval':critical_value, 'pval': pval}

    return result
for measure in aan_global_measures.columns: 
    eff = monte_carlo_2(aan_global_measures[measure].values, hc_global_measures[measure].values, permutations=50000)
    print(f'Group difference between AAN and HC for {measure}','pval=', eff['pval'],'critical value:',eff['critval'])
for measure in aan_global_measures.columns: 
    eff = monte_carlo_2(aan_global_measures[measure].values, wr_global_measures[measure].values, permutations=50000)
    print(f'Group difference between AAN and WR for {measure}','pval=', eff['pval'],'critical value:',eff['critval'])
for measure in aan_global_measures.columns: 
    eff = monte_carlo_2(hc_global_measures[measure].values, wr_global_measures[measure].values, permutations=50000)
    print(f'Group difference between HC and WR for {measure}','pval=', eff['pval'],'critical value:',eff['critval'])<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1310016738613338, 'fraction_numerical': 0.011012245617126245, 'mean_word_length': 4.03861517976032, 'pattern_counts': {'":': 0, '<': 14, '<?xml version=': 0, '>': 22, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '2', 'dir': 'jupyter-structured-clean-dedup', 'id': '609927', 'license': 'MIT', 'n_tokens_mistral': 3933, 'n_tokens_neox': 3667, 'n_words': 787, 'path': 'scripts/work_in_progress_scripts/group_differences_on_global_measures.ipynb', 'repo_name': 'WMDA/SCN', 'size': '11351'}