import numpy as np
import pandas as pd
import scipy.stats as stats
import matplotlib.pyplot as plt
import seaborn as sns

# 设置中文显示
plt.rcParams['font.sans-serif'] = ['SimHei']  # 用来正常显示中文标签
plt.rcParams['axes.unicode_minus'] = False    # 用来正常显示负号

# 数据输入
A1 = [0.88, 0.85, 0.79, 0.86, 0.85, 0.83]
A2 = [0.87, 0.92, 0.85, 0.83, 0.90, 0.80]
A3 = [0.84, 0.78, 0.81, 0.80, 0.85, 0.83]
A4 = [0.81, 0.86, 0.90, 0.87, 0.78, 0.79]

# 创建数据框
data = pd.DataFrame({
    'A1': A1,
    'A2': A2,
    'A3': A3,
    'A4': A4
})

# 计算每种催化剂的均值和标准差
means = [np.mean(A1), np.mean(A2), np.mean(A3), np.mean(A4)]
stds = [np.std(A1, ddof=1), np.std(A2, ddof=1), np.std(A3, ddof=1), np.std(A4, ddof=1)]

print("各催化剂的得率均值:")
for i, mean in enumerate(means):
    print(f"A{i+1}: {mean:.4f}")

print("\n各催化剂的得率标准差:")
for i, std in enumerate(stds):
    print(f"A{i+1}: {std:.4f}")

# 进行单因素方差分析
# 将数据转换为适合ANOVA的格式
df_melt = pd.melt(data)
df_melt.columns = ['催化剂', '得率']

# 执行ANOVA
f_stat, p_value = stats.f_oneway(A1, A2, A3, A4)

print(f"\n方差分析结果:")
print(f"F统计量: {f_stat:.4f}")
print(f"P值: {p_value:.4f}")

alpha = 0.05
if p_value < alpha:
    conclusion = "拒绝原假设，四种催化剂对产品得率有显著影响"
else:
    conclusion = "接受原假设，四种催化剂对产品得率无显著影响"

print(f"\n在显著性水平α={alpha}下，{conclusion}")

# 绘制箱线图
plt.figure(figsize=(10, 6))
sns.boxplot(x='催化剂', y='得率', data=df_melt)
plt.title('四种催化剂对产品得率的影响', fontsize=14)
plt.grid(True, linestyle='--', alpha=0.7)
plt.ylim(0.75, 0.95)  # 调整y轴范围使差异更明显
plt.tight_layout()
plt.show()

# 如果ANOVA结果显著，进行事后多重比较（图基检验）
if p_value < alpha:
    from statsmodels.stats.multicomp import pairwise_tukeyhsd
    print("\n事后多重比较(Tukey HSD)结果:")
    tukey = pairwise_tukeyhsd(df_melt['得率'], df_melt['催化剂'], alpha=alpha)
    print(tukey)

# 进一步检验ANOVA的假设条件
# 1. 正态性检验
print("\n正态性检验(Shapiro-Wilk)结果:")
for i, group in enumerate([A1, A2, A3, A4]):
    stat, p = stats.shapiro(group)
    print(f"A{i+1}: 统计量={stat:.4f}, P值={p:.4f}, {'正态分布' if p > 0.05 else '非正态分布'}")

# 2. 方差齐性检验
stat, p = stats.levene(A1, A2, A3, A4)
print(f"\n方差齐性检验(Levene)结果: 统计量={stat:.4f}, P值={p:.4f}")
print(f"{'方差齐性' if p > 0.05 else '方差不齐'}")

# 描述统计量表格
stats_df = pd.DataFrame({
    '催化剂': ['A1', 'A2', 'A3', 'A4'],
    '样本量': [len(A1), len(A2), len(A3), len(A4)],
    '均值': means,
    '标准差': stds,
    '最小值': [min(A1), min(A2), min(A3), min(A4)],
    '最大值': [max(A1), max(A2), max(A3), max(A4)]
})

print("\n描述统计量:")
print(stats_df.to_string(index=False))