# import naginterfaces.library.anova as aov
import statsmodels.stats.multicomp as mc
import scipy.stats as st
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

dfl = pd.DataFrame({'feed':['A', 'B', 'C', 'D']*4,
          'crop': [60, 62, 68, 56, 59, 67, 67, 60,
                    59, 63, 71, 59, 58, 64, 70, 57]})

comp = mc.MultiComparison(dfl['crop'], dfl['feed'])

hsd = comp.tukeyhsd()
print('Tukey HSD:\n', hsd.summary())

bonf, a1, b1 = comp.allpairtest(st.ttest_ind, method= "bonf")
print('Bonferoni:\n', bonf)

sidak, a2, b2 = comp.allpairtest(st.ttest_ind, method= "sidak")
print('Sidak:\n', sidak)

scheffe, a3, b3 = comp.allpairtest(st.ttest_ind, method= "sh")
print('Scheffe:\n', scheffe)

#---------- HSD 그림 ------------------------------
fig, ax1 = plt.subplots(1,1)
hsd.plot_simultaneous(ax=ax1)
mm = np.mean(dfl['crop'])
ax1.axvline(x=mm, color='r', label='gross mean at: {}'.format(mm))
ax1.legend()
######## 등분산 독립 2표본 검정결과 비교 #########################
dfw = pd.DataFrame({'A':[60,59,59,58], 'B':[62,67,63,64],
                    'C':[68,67,71,70], 'D':[56,60,59,57]})

t1, p1 = st.ttest_ind(dfw['A'], dfw['B'])
t2, p2 = st.ttest_ind(dfw['A'], dfw['C'])
print(t1, p1)
print(t2, p2)