Updated Statistics file conducting cross labeling statisitcal analysis

import matplotlib.pyplot as plt

import seaborn as sns

from scipy.stats import chi2_contingency, shapiro, levene, ttest_ind, mannwhitneyu

from statsmodels.stats.multitest import multipletests

# Load Excel file

custom_colors = ["#38363A", "#BAAC97", "#902C44"]

def load_data(file_path, sheet_name):

    return pd.read_excel(file_path, sheet_name=sheet_name)

# Data Preparation

def preprocess_data(data):

    data = data[data['GL Optimized CV'].isin([1, 2])]

    data['final_DX'] = data['final_DX'].replace({

        1: "CN", 2: "Dementia due to AD", 3: "MCI due to AD",

        5: "MCI due to other etiology", -8: "MCI due to AD"

    })

    final_dx_labels = {1: "CN", 2: "Dementia due to AD", 3: "MCI due to AD", 5: "MCI due to other etiology",

                       -8: "MCI due to AD"}

    data['final_DX'] = data['final_DX'].map(final_dx_labels)

    data['AAB'] = pd.to_numeric(data['AAB'], errors='coerce')

    return data

# Crosstab with percentages

def calculate_percentage_table(data, row_var, col_var):

    return pd.crosstab(data[row_var], data[col_var], normalize='index') * 100

# Perform Chi-squared tests

def chi_square_test(data, col):

    table = pd.crosstab(data[col], data['GL Optimized CV'])

    chi2, p, _, _ = chi2_contingency(table)

    return p

# Summary statistics for Age at Baseline

def summarize_column(data, column):

    return data[column].describe()

# Convert variables to numeric

def convert_to_numeric(data, variables):

def perform_chi_square_tests(data, variables):

    for var in variables:

        if var in data.columns:

            data[var] = pd.to_numeric(data[var], errors='coerce')

    return data

        contingency_table = pd.crosstab(data[var], data['GL Optimized CV'])

        chi2, p, _, _ = chi2_contingency(contingency_table)

        print(f"Chi-square test for {var}: chi2={chi2}, p-value={p}")

# Statistical tests (T-test or Mann-Whitney U test)

def perform_statistical_tests(data, variables):

def calculate_statistics(data, variables):

    results = []

    for var in variables:

        for g1, g2 in [(1, 2)]:

        if var in data.columns:

            for g1 in range(1, 3):

                for g2 in range(g1 + 1, 4):

                    group1 = data[data['GL Optimized CV'] == g1][var].dropna()

                    group2 = data[data['GL Optimized CV'] == g2][var].dropna()

                    if len(group1) < 15 or len(group2) < 15:

                u_stat, p_val = mannwhitneyu(group1, group2, alternative='less')

                test_type = "Mann-Whitney U test"

                        stat, p_value = mannwhitneyu(group1, group2, alternative='two-sided')

                        test_type = "Mann-Whitney U"

                    else:

                norm_g1 = shapiro(group1).pvalue > 0.05

                norm_g2 = shapiro(group2).pvalue > 0.05

                        normal_g1 = shapiro(group1)[1] > 0.05

                        normal_g2 = shapiro(group2)[1] > 0.05

                        equal_var = levene(group1, group2)[1] > 0.05

                levene_p = levene(group1, group2).pvalue

                equal_var = levene_p > 0.05

                if norm_g1 and norm_g2:

                    t_stat, p_val = ttest_ind(group1, group2, equal_var=equal_var, alternative="less")

                        if normal_g1 and normal_g2:

                            stat, p_value = ttest_ind(group1, group2, equal_var=equal_var)

                            test_type = "Student's t-test" if equal_var else "Welch's t-test"

                        else:

                    t_stat, p_val = ttest_ind(group1, group2, equal_var=equal_var, alternative="less")

                    test_type = "Non-parametric t-test (due to normality issues)"

                            stat, p_value = mannwhitneyu(group1, group2, alternative='two-sided')

                            test_type = "Mann-Whitney U"

            results.append([var, f"Group {g1} vs {g2}", test_type, p_val])

                    results.append([var, f"Group {g1} vs {g2}", test_type, round(stat, 3), round(p_value, 4)])

    results_df = pd.DataFrame(results, columns=["Variable", "Comparison", "Test_Type", "p_value"])

    results_df['corrected_p_value'] = multipletests(results_df['p_value'], method='bonferroni')[1]

    return results_df

    results_df = pd.DataFrame(results, columns=["Variable", "Comparison", "Test Type", "Statistic", "p-value"])

    results_df["Corrected p-value"] = np.minimum(results_df["p-value"] * len(results_df), 1.0)

    print(results_df)

# Plotting boxplots

def plot_boxplots(data, variables, title_map):

    for var in variables:

def boxplot_variable(data, var, ylabel):

    plt.figure(figsize=(8, 5))

        sns.boxplot(x='GL Optimized CV', y=var, hue='final_DX', data=data, palette=["#38363A", "#BAAC97", "#902C44"])

    sns.boxplot(x='GL Optimized CV', y=var, hue='final_DX', data=data, palette=custom_colors)

    plt.xlabel("GL Optimized CV")

        plt.ylabel(title_map.get(var, var))

    plt.ylabel(ylabel)

    plt.title(f"Boxplot for {ylabel} by community")

    plt.legend(title="Diagnosis")

        plt.grid(False)

    # plt.savefig(f"boxplot_{var}.png", dpi=300)

    plt.show()

# Main execution

if __name__ == "__main__":

    file_path = "/Users/samansarraf/Desktop/Stanford/Data_for_statistics.xlsx"

    sheet_name = "SheetName"

    data = load_data(file_path, sheet_name)

    data = preprocess_data(data)

    print("Diagnosis Table:")

    print(pd.crosstab(data['final_DX'], data['GL Optimized CV']))

    print("\nPercentage Table:")

    print(calculate_percentage_table(data, 'final_DX', 'GL Optimized CV'))

    print("\nSex Count:")

    print(pd.crosstab(data['final_DX'], data['GL Optimized CV'], normalize='index') * 100)

    print(data['SEX'].value_counts())

    for col in ['APOE', 'SEX']:

        print(f"\n{col} Chi-squared test p-value: {chi_square_test(data, col):.4f}")

    print("\nAge at Baseline Summary:")

    print(summarize_column(data, 'AAB'))

    variables = ["AAB", "EDUC", "SEX", "MMSE", "MOCA", "ADAS13", "CDRSB",

                 "AV45_cs", "ABETA", "TAU", "PTAU", "Hippocampus", "Entorhinal"]

    data = convert_to_numeric(data, variables)

    results_df = perform_statistical_tests(data, variables)

    print("\nStatistical Test Results:")

    print(results_df)

    print("\nAmyloid Positivity Table:")

    print(pd.crosstab(data['AV45_cs_pos'], data['GL Optimized CV']))

    perform_chi_square_tests(data, ['APOE', 'SEX'])

    plot_variables = ["AV45_cs", "TAU", "PTAU", "Hippocampus"]

    title_map = {"AV45_cs": "Aβ PET", "TAU": "CSF TAU", "PTAU": "CSF PTAU", "Hippocampus": "Hippocampal Volume"}

    plot_boxplots(data, plot_variables, title_map)

    variables = ["AAB", "SEX", "MMSE", "MOCA", "ADAS13", "CDRSB", "AV45_cs", "ABETA", "TAU", "PTAU", "Hippocampus",

                 "Entorhinal"]

    calculate_statistics(data, variables)

    for var in ["MMSE", "MOCA", "CDRSB", "ADAS13", "PTAU", "TAU", "AV45_cs", "Hippocampus"]:

        boxplot_variable(data, var, var)