added Statistics file conducting cross labeling statisitcal analysis

import pandas as pd

import numpy as np

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

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"

    })

    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):

    for var in variables:

        if var in data.columns:

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

    return data

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

def perform_statistical_tests(data, variables):

    results = []

    for var in variables:

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

            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"

            else:

                norm_g1 = shapiro(group1).pvalue > 0.05

                norm_g2 = shapiro(group2).pvalue > 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")

                    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)"

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

    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

# Plotting boxplots

def plot_boxplots(data, variables, title_map):

    for var in variables:

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

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

        plt.xlabel("GL Optimized CV")

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

        plt.legend(title="Diagnosis")

        plt.grid(False)

        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(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']))

    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)