<html><head><meta name="color-scheme" content="light dark"></head><body><pre style="word-wrap: break-word; white-space: pre-wrap;">import numpy as np
import os
import matplotlib.pyplot as plt
import scipy
from scipy import optimize
from scipy import stats
import math
import pandas
import skimage.measure as sk
from statsmodels.stats.multicomp import pairwise_tukeyhsd

def main():

    directoryname = '/home/morrile2/Documents/Multis/app/UltrasoundThickness/TrainingData/'

    F_M_bound = []
    M_B_bound = []
    S_F_bound = []
    T_S_bound = []

    for dirName, subdirList, fileList in os.walk(directoryname):
        for filename in fileList:
            if "_fm.npy" in filename.lower():  # check whether the file's Fat/Muscle boundary
                F_M_bound.append(os.path.join(dirName, filename))
            elif "_mb.npy" in filename.lower():  # check whether the file's Muscle/Bone boundary
                M_B_bound.append(os.path.join(dirName, filename))
            elif "_sf.npy" in filename.lower():  # check whether the file's Skin/Fat boundary
                S_F_bound.append(os.path.join(dirName, filename))
            elif "_ts.npy" in filename.lower():  # check whether the file's Transducer/Skin boundary
                T_S_bound.append(os.path.join(dirName, filename))

    TS_res = []
    SF_res = []
    FM_res = []
    MB_res = []

    TS_inv = []
    SF_inv = []
    FM_inv = []
    MB_inv = []


    pp = 0

    print("T/S")
    figure = plt.figure(1)
    figure.add_subplot(221)
    for i in T_S_bound:
        ydist, length = getYDist(np.load(i))
        x = np.linspace(-3,3,length)
        TS_res.append(curveFit(x, ydist))
        # TS_res.append(getYDist(np.load(i)))
        m = sk.moments(np.load(i)/np.mean(np.load(i)))
        m_c = sk.moments_central(np.load(i)/np.mean(np.load(i)), m[0,1]/m[0,0], m[1,0]/m[0,0])
        TS_inv.append(sk.moments_hu(m_c)[pp])
        plt.plot(x, ydist)
    plt.title("Transducer/Skin")

    print("S/F")
    figure.add_subplot(222)
    for i in S_F_bound:
        ydist, length = getYDist(np.load(i))
        x = np.linspace(-3,3,length)
        SF_res.append(curveFit(x, ydist))
        # SF_res.append(getYDist(np.load(i)))
        m = sk.moments(np.load(i)/np.mean(np.load(i)))
        m_c = sk.moments_central(np.load(i)/np.mean(np.load(i)), m[0,1]/m[0,0], m[1,0]/m[0,0])
        SF_inv.append(sk.moments_hu(m_c)[pp])
        plt.plot(x, ydist)
    plt.title("Skin/Fat")

    print("F/M")
    figure.add_subplot(223)
    for i in F_M_bound:
        ydist, length = getYDist(np.load(i))
        x = np.linspace(-3,3,length)
        FM_res.append(curveFit(x, ydist))
        # FM_res.append(getYDist(np.load(i)))
        m = sk.moments(np.load(i)/np.mean(np.load(i)))
        m_c = sk.moments_central(np.load(i)/np.mean(np.load(i)), m[0,1]/m[0,0], m[1,0]/m[0,0])
        FM_inv.append(sk.moments_hu(m_c)[pp])
        plt.plot(x, ydist)
    plt.title("Fat/Muscle")

    print("M/B")
    figure.add_subplot(224)
    for i in M_B_bound:
        ydist, length = getYDist(np.load(i))
        x = np.linspace(-3,3,length)
        MB_res.append(curveFit(x, ydist))
        # MB_res.append(getYDist(np.load(i)))
        m = sk.moments(np.load(i)/np.mean(np.load(i)))
        m_c = sk.moments_central(np.load(i)/np.mean(np.load(i)), m[0,1]/m[0,0], m[1,0]/m[0,0])
        MB_inv.append(sk.moments_hu(m_c)[pp])
        plt.plot(x, ydist)
    plt.legend
    plt.title("Muscle/Bone")

    #Create dataframe of mean and standard deviations
    df = pandas.DataFrame(data=None, columns=['TS_mean', 'TS_std', 'SF_mean', 'SF_std', 'FM_mean', 'FM_std', 'MB_mean', 'MB_std', 'TS_I1', 'SF_I1', 'FM_I1', 'MB_I1'])
    df_TS = pandas.DataFrame(np.array(zip(*TS_res)).T, columns=['TS_mean', 'TS_std'])
    df_SF = pandas.DataFrame(np.array(zip(*SF_res)).T, columns=['SF_mean', 'SF_std'])
    df_FM = pandas.DataFrame(np.array(zip(*FM_res)).T, columns=['FM_mean', 'FM_std'])
    df_MB = pandas.DataFrame(np.array(zip(*MB_res)).T, columns=['MB_mean', 'MB_std'])

    df_TSmom = pandas.DataFrame(np.array(TS_inv).T, columns=['TS_I1'])
    df_SFmom = pandas.DataFrame(np.array(SF_inv).T, columns=['SF_I1'])
    df_FMmom = pandas.DataFrame(np.array(FM_inv).T, columns=['FM_I1'])
    df_MBmom = pandas.DataFrame(np.array(MB_inv).T, columns=['MB_I1'])

    df = df.append(pandas.concat([df_TS, df_SF, df_FM, df_MB, df_TSmom, df_SFmom, df_FMmom, df_MBmom], axis=1), ignore_index=True)

    print(df.mean())
    f2 = plt.figure(2)
    f2.add_subplot(121)
    df.boxplot(column = ['TS_mean', 'SF_mean', 'FM_mean', 'MB_mean'], return_type='dict')
    f2.add_subplot(122)
    df.boxplot(column = ['TS_std', 'SF_std', 'FM_std', 'MB_std'], return_type='dict')

    print("Mean")
    f_val, p_val = stats.f_oneway(df['TS_mean'], df['SF_mean'], df['FM_mean'], df['MB_mean'])
    print(f_val, p_val)

    groups = []
    data = []

    for k in df['TS_mean']:
        data.append(k)
        groups.append('TS')
    for k in df['SF_mean']:
        data.append(k)
        groups.append('SF')
    for k in df['FM_mean']:
        data.append(k)
        groups.append('FM')
    for k in df['MB_mean']:
        data.append(k)
        groups.append('MB')

    tukey = pairwise_tukeyhsd(endog=np.array(data), groups = np.array(groups), alpha=0.05)
    print(tukey)

    groups = []
    data = []
    print("Standard Deviation")
    for k in df['TS_std']:
        data.append(k)
        groups.append('TS')
    for k in df['SF_std']:
        data.append(k)
        groups.append('SF')
    for k in df['FM_std']:
        data.append(k)
        groups.append('FM')
    for k in df['MB_std']:
        data.append(k)
        groups.append('MB')

    tukey = pairwise_tukeyhsd(endog=np.array(data), groups = np.array(groups), alpha=0.05)
    print(tukey)


    plt.figure(3)
    df.boxplot(column=['TS_I1', 'SF_I1', 'FM_I1', 'MB_I1'], return_type='dict')


    f_val, p_val = stats.f_oneway(df['TS_std'], df['SF_std'], df['FM_std'], df['MB_std'])
    print(f_val, p_val)

    print('First Invariant')
    f_val, p_val = stats.f_oneway(df['TS_I1'], df['SF_I1'], df['FM_I1'], df['MB_I1'])
    print(f_val, p_val)

    plt.show()

def getYDist(Data_mat):
    # p = Data_mat.shape[0]
    # mean = np.mean(Data_mat[int(p/4):int(3*p/4), int(p/8):int(3*p/8)], axis=1)
    mean = np.mean(Data_mat, axis=1)
    mean = mean/np.max(mean)
    return mean, len(mean)
    # return(np.mean(np.gradient(mean)), np.cov(mean))

def curveFit(x, y):
    popt, pcov = scipy.optimize.curve_fit(gaussCurve, x, y, p0=[0, .1])
    return popt

def gaussCurve(x, mu, sigma):
    return 1/(sigma*math.sqrt(2*math.pi)) * np.exp(-(x - mu) ** 2 / (2 * sigma ** 2))

if __name__ == "__main__":
    main()</pre></body></html>