Python实例题:Python计算数理统计
目录
Python实例题
题目
代码实现
实现原理
描述性统计:
参数估计:
假设检验:
相关与回归分析:
可视化功能:
关键代码解析
1. 描述性统计
2. 置信区间计算
3. 假设检验
4. 线性回归分析
使用说明
安装依赖:
基本用法:
示例输出:
扩展建议
增强功能:
用户界面:
性能优化:
教学辅助:
Python实例题
题目
Python计算数理统计
代码实现
import numpy as np
import matplotlib.pyplot as plt
from scipy import stats
import seaborn as sns
import pandas as pd
from sympy import symbols, simplify, latex
from statsmodels.formula.api import ols
import statsmodels.api as smclass StatisticalAnalysis:"""数理统计分析类,支持描述性统计、参数估计和假设检验"""def __init__(self):"""初始化分析器"""passdef descriptive_statistics(self, data):"""计算描述性统计量参数:data: 样本数据(列表或NumPy数组)返回:dict: 包含各种描述性统计量的字典"""results = {'样本大小': len(data),'均值': np.mean(data),'中位数': np.median(data),'众数': stats.mode(data, keepdims=True)[0][0],'标准差': np.std(data, ddof=1), # 无偏估计'方差': np.var(data, ddof=1), # 无偏估计'最小值': np.min(data),'最大值': np.max(data),'极差': np.max(data) - np.min(data),'四分位数Q1': np.percentile(data, 25),'四分位数Q2': np.percentile(data, 50),'四分位数Q3': np.percentile(data, 75),'四分位距': np.percentile(data, 75) - np.percentile(data, 25),'偏度': stats.skew(data),'峰度': stats.kurtosis(data, fisher=True) # Fisher校正的峰度}return resultsdef confidence_interval(self, data, confidence_level=0.95):"""计算均值的置信区间参数:data: 样本数据(列表或NumPy数组)confidence_level: 置信水平,默认为0.95返回:tuple: (置信下限, 置信上限)"""n = len(data)mean = np.mean(data)std = np.std(data, ddof=1) # 样本标准差alpha = 1 - confidence_level# 计算t分布的临界值t_critical = stats.t.ppf(1 - alpha/2, df=n-1)# 计算标准误差standard_error = std / np.sqrt(n)# 计算置信区间lower = mean - t_critical * standard_errorupper = mean + t_critical * standard_errorreturn (lower, upper)def hypothesis_test_one_sample(self, data, mu0, alternative='two-sided', alpha=0.05):"""单样本t检验参数:data: 样本数据(列表或NumPy数组)mu0: 原假设的总体均值alternative: 备择假设方向,'two-sided'(双侧)、'greater'(大于)或'less'(小于)alpha: 显著性水平,默认为0.05返回:dict: 包含检验统计量、p值和检验结论的字典"""# 执行t检验t_stat, p_value = stats.ttest_1samp(data, mu0)# 根据备择假设调整p值if alternative == 'greater' and t_stat > 0:p_value = p_value / 2elif alternative == 'less' and t_stat < 0:p_value = p_value / 2elif alternative in ['greater', 'less']:p_value = 1 - p_value / 2# 得出结论reject_null = p_value < alphaconclusion = f"在显著性水平 {alpha} 下,{'拒绝' if reject_null else '不拒绝'} 原假设"results = {'检验统计量': t_stat,'p值': p_value,'自由度': len(data) - 1,'结论': conclusion}return resultsdef hypothesis_test_two_samples(self, data1, data2, equal_var=True, alternative='two-sided', alpha=0.05):"""两样本t检验参数:data1: 第一个样本数据data2: 第二个样本数据equal_var: 是否假设两总体方差相等,默认为Truealternative: 备择假设方向,'two-sided'(双侧)、'greater'(大于)或'less'(小于)alpha: 显著性水平,默认为0.05返回:dict: 包含检验统计量、p值和检验结论的字典"""# 执行t检验t_stat, p_value = stats.ttest_ind(data1, data2, equal_var=equal_var)# 根据备择假设调整p值if alternative == 'greater' and t_stat > 0:p_value = p_value / 2elif alternative == 'less' and t_stat < 0:p_value = p_value / 2elif alternative in ['greater', 'less']:p_value = 1 - p_value / 2# 得出结论reject_null = p_value < alphaconclusion = f"在显著性水平 {alpha} 下,{'拒绝' if reject_null else '不拒绝'} 原假设"results = {'检验统计量': t_stat,'p值': p_value,'自由度': len(data1) + len(data2) - 2 if equal_var else None,'结论': conclusion}return resultsdef anova_one_way(self, *groups):"""单因素方差分析参数:*groups: 多个样本组(每个组是一个列表或NumPy数组)返回:dict: 包含F统计量、p值和检验结论的字典"""# 执行单因素方差分析f_stat, p_value = stats.f_oneway(*groups)# 得出结论alpha = 0.05reject_null = p_value < alphaconclusion = f"在显著性水平 {alpha} 下,{'至少有一组均值存在显著差异' if reject_null else '各组均值无显著差异'}"results = {'F统计量': f_stat,'p值': p_value,'结论': conclusion}return resultsdef correlation_analysis(self, X, Y, method='pearson'):"""相关性分析参数:X: 第一个变量的数据Y: 第二个变量的数据method: 相关系数计算方法,'pearson'(皮尔逊)或'spearman'(斯皮尔曼)返回:dict: 包含相关系数、p值和结论的字典"""if method == 'pearson':corr, p_value = stats.pearsonr(X, Y)elif method == 'spearman':corr, p_value = stats.spearmanr(X, Y)else:raise ValueError("不支持的相关系数计算方法")# 得出结论alpha = 0.05significant = p_value < alphastrength = ''if abs(corr) < 0.2:strength = '极弱相关'elif abs(corr) < 0.4:strength = '弱相关'elif abs(corr) < 0.6:strength = '中等相关'elif abs(corr) < 0.8:strength = '强相关'else:strength = '极强相关'direction = '正相关' if corr > 0 else '负相关'results = {f'{method}相关系数': corr,'p值': p_value,'显著性': f"在显著性水平 {alpha} 下,{'显著' if significant else '不显著'}",'相关强度和方向': f"{strength},{direction}"}return resultsdef linear_regression(self, X, Y):"""简单线性回归分析参数:X: 自变量数据Y: 因变量数据返回:dict: 包含回归系数、截距、R²值和p值的字典"""# 添加常数项X_with_const = sm.add_constant(X)# 执行线性回归model = sm.OLS(Y, X_with_const).fit()# 提取结果results = {'截距': model.params[0],'斜率': model.params[1],'R²值': model.rsquared,'调整R²值': model.rsquared_adj,'F统计量': model.fvalue,'F检验p值': model.f_pvalue,'斜率p值': model.pvalues[1],'回归方程': f"Y = {model.params[0]:.4f} + {model.params[1]:.4f}X",'模型摘要': model.summary()}return resultsdef chi_square_test(self, observed):"""卡方检验(拟合优度检验或独立性检验)参数:observed: 观察频数的二维数组或列表的列表返回:dict: 包含卡方统计量、p值、自由度和结论的字典"""# 执行卡方检验chi2, p, dof, expected = stats.chi2_contingency(observed)# 得出结论alpha = 0.05reject_null = p < alphaconclusion = f"在显著性水平 {alpha} 下,{'拒绝' if reject_null else '不拒绝'} 原假设"results = {'卡方统计量': chi2,'p值': p,'自由度': dof,'期望频数': expected,'结论': conclusion}return resultsdef plot_histogram(self, data, bins=30, title=None):"""绘制样本直方图参数:data: 样本数据bins: 直方图的箱数,默认为30title: 图像标题,默认为None"""plt.figure(figsize=(10, 6))sns.histplot(data, bins=bins, kde=True)plt.grid(True)plt.axvline(x=np.mean(data), color='r', linestyle='--', label=f'均值: {np.mean(data):.2f}')plt.axvline(x=np.median(data), color='g', linestyle='--', label=f'中位数: {np.median(data):.2f}')if title:plt.title(title, fontsize=14)else:plt.title(f"样本直方图", fontsize=14)plt.xlabel('值', fontsize=12)plt.ylabel('频率', fontsize=12)plt.legend()plt.show()def plot_boxplot(self, data, title=None):"""绘制箱线图参数:data: 样本数据(可以是单个列表或多个列表的列表)title: 图像标题,默认为None"""plt.figure(figsize=(10, 6))if isinstance(data[0], list) or isinstance(data[0], np.ndarray):# 多个数据集plt.boxplot(data)plt.xticks(range(1, len(data) + 1), [f"数据{i+1}" for i in range(len(data))])else:# 单个数据集plt.boxplot(data)plt.grid(True)if title:plt.title(title, fontsize=14)else:plt.title(f"箱线图", fontsize=14)plt.ylabel('值', fontsize=12)plt.show()def plot_scatter(self, X, Y, title=None, xlabel=None, ylabel=None):"""绘制散点图参数:X: 自变量数据Y: 因变量数据title: 图像标题,默认为Nonexlabel: x轴标签,默认为Noneylabel: y轴标签,默认为None"""plt.figure(figsize=(10, 6))plt.scatter(X, Y, alpha=0.7)# 计算并绘制回归线slope, intercept, r_value, p_value, std_err = stats.linregress(X, Y)plt.plot(X, intercept + slope*X, 'r-', label=f'回归线: y = {intercept:.4f} + {slope:.4f}x')plt.grid(True)if title:plt.title(title, fontsize=14)else:plt.title(f"散点图 (r={r_value:.4f})", fontsize=14)if xlabel:plt.xlabel(xlabel, fontsize=12)else:plt.xlabel('X', fontsize=12)if ylabel:plt.ylabel(ylabel, fontsize=12)else:plt.ylabel('Y', fontsize=12)plt.legend()plt.show()# 示例使用
def example_usage():analysis = StatisticalAnalysis()print("\n===== 描述性统计示例 =====")# 生成样本数据np.random.seed(42) # 设置随机种子,确保结果可重现data = np.random.normal(loc=50, scale=10, size=100)# 计算描述性统计量stats = analysis.descriptive_statistics(data)for key, value in stats.items():print(f"{key}: {value:.4f}" if isinstance(value, float) else f"{key}: {value}")# 绘制直方图analysis.plot_histogram(data, title="样本数据直方图")# 绘制箱线图analysis.plot_boxplot(data, title="样本数据箱线图")print("\n===== 参数估计示例 =====")# 计算置信区间lower, upper = analysis.confidence_interval(data, confidence_level=0.95)print(f"95%置信区间: [{lower:.4f}, {upper:.4f}]")print("\n===== 假设检验示例 =====")# 单样本t检验test_results = analysis.hypothesis_test_one_sample(data, mu0=50)for key, value in test_results.items():print(f"{key}: {value}")# 两样本t检验data2 = np.random.normal(loc=55, scale=10, size=100)test_results = analysis.hypothesis_test_two_samples(data, data2, equal_var=True)for key, value in test_results.items():print(f"{key}: {value}")# 单因素方差分析group1 = np.random.normal(loc=50, scale=10, size=50)group2 = np.random.normal(loc=55, scale=10, size=50)group3 = np.random.normal(loc=60, scale=10, size=50)anova_results = analysis.anova_one_way(group1, group2, group3)for key, value in anova_results.items():print(f"{key}: {value}")print("\n===== 相关性分析示例 =====")# 生成相关数据X = np.random.normal(0, 1, 100)Y = 2*X + np.random.normal(0, 1, 100)# 计算相关系数corr_results = analysis.correlation_analysis(X, Y, method='pearson')for key, value in corr_results.items():print(f"{key}: {value}")# 绘制散点图analysis.plot_scatter(X, Y, title="X和Y的散点图", xlabel="X", ylabel="Y")print("\n===== 线性回归分析示例 =====")# 执行线性回归regression_results = analysis.linear_regression(X, Y)for key, value in regression_results.items():if key != '模型摘要':print(f"{key}: {value}")print("\n===== 卡方检验示例 =====")# 创建列联表observed = np.array([[20, 30, 10], [15, 25, 20]])chi2_results = analysis.chi_square_test(observed)for key, value in chi2_results.items():if key != '期望频数':print(f"{key}: {value}")else:print(f"{key}:\n{value}")if __name__ == "__main__":example_usage()
实现原理
这个数理统计计算工具基于以下技术实现:
-
描述性统计:
- 计算均值、中位数、众数、标准差等基本统计量
- 分析数据分布特征(偏度、峰度)
- 提供直观的数据可视化(直方图、箱线图)
-
参数估计:
- 计算样本均值的置信区间
- 基于 t 分布进行区间估计
- 支持不同置信水平的选择
-
假设检验:
- 实现单样本和两样本 t 检验
- 支持单因素方差分析(ANOVA)
- 提供卡方检验(拟合优度和独立性检验)
-
相关与回归分析:
- 计算皮尔逊和斯皮尔曼相关系数
- 执行简单线性回归分析
- 评估回归模型的拟合优度和显著性
-
可视化功能:
- 使用 Matplotlib 和 Seaborn 绘制统计图表
- 直观展示数据分布和变量关系
- 支持自定义图表样式和标题
关键代码解析
1. 描述性统计
def descriptive_statistics(self, data):"""计算描述性统计量"""results = {'样本大小': len(data),'均值': np.mean(data),'中位数': np.median(data),'众数': stats.mode(data, keepdims=True)[0][0],'标准差': np.std(data, ddof=1),'方差': np.var(data, ddof=1),'最小值': np.min(data),'最大值': np.max(data),'极差': np.max(data) - np.min(data),'四分位数Q1': np.percentile(data, 25),'四分位数Q2': np.percentile(data, 50),'四分位数Q3': np.percentile(data, 75),'四分位距': np.percentile(data, 75) - np.percentile(data, 25),'偏度': stats.skew(data),'峰度': stats.kurtosis(data, fisher=True)}return results
2. 置信区间计算
def confidence_interval(self, data, confidence_level=0.95):"""计算均值的置信区间"""n = len(data)mean = np.mean(data)std = np.std(data, ddof=1)alpha = 1 - confidence_levelt_critical = stats.t.ppf(1 - alpha/2, df=n-1)standard_error = std / np.sqrt(n)lower = mean - t_critical * standard_errorupper = mean + t_critical * standard_errorreturn (lower, upper)
3. 假设检验
def hypothesis_test_one_sample(self, data, mu0, alternative='two-sided', alpha=0.05):"""单样本t检验"""t_stat, p_value = stats.ttest_1samp(data, mu0)if alternative == 'greater' and t_stat > 0:p_value = p_value / 2elif alternative == 'less' and t_stat < 0:p_value = p_value / 2elif alternative in ['greater', 'less']:p_value = 1 - p_value / 2reject_null = p_value < alphaconclusion = f"在显著性水平 {alpha} 下,{'拒绝' if reject_null else '不拒绝'} 原假设"results = {'检验统计量': t_stat,'p值': p_value,'自由度': len(data) - 1,'结论': conclusion}return results
4. 线性回归分析
def linear_regression(self, X, Y):"""简单线性回归分析"""X_with_const = sm.add_constant(X)model = sm.OLS(Y, X_with_const).fit()results = {'截距': model.params[0],'斜率': model.params[1],'R²值': model.rsquared,'调整R²值': model.rsquared_adj,'F统计量': model.fvalue,'F检验p值': model.f_pvalue,'斜率p值': model.pvalues[1],'回归方程': f"Y = {model.params[0]:.4f} + {model.params[1]:.4f}X",'模型摘要': model.summary()}return results
使用说明
-
安装依赖:
pip install numpy matplotlib scipy seaborn pandas statsmodels
-
基本用法:
from statistical_analysis import StatisticalAnalysis# 创建分析器实例
analysis = StatisticalAnalysis()# 生成样本数据
data = np.random.normal(loc=50, scale=10, size=100)# 计算描述性统计量
stats = analysis.descriptive_statistics(data)
for key, value in stats.items():print(f"{key}: {value}")# 计算置信区间
lower, upper = analysis.confidence_interval(data, confidence_level=0.95)
print(f"95%置信区间: [{lower:.4f}, {upper:.4f}]")# 单样本t检验
test_results = analysis.hypothesis_test_one_sample(data, mu0=50)
for key, value in test_results.items():print(f"{key}: {value}")# 相关性分析
X = np.random.normal(0, 1, 100)
Y = 2*X + np.random.normal(0, 1, 100)
corr_results = analysis.correlation_analysis(X, Y)
for key, value in corr_results.items():print(f"{key}: {value}")# 线性回归分析
regression_results = analysis.linear_regression(X, Y)
for key, value in regression_results.items():if key != '模型摘要':print(f"{key}: {value}")# 绘制散点图
analysis.plot_scatter(X, Y, title="X和Y的散点图")
-
示例输出:
样本大小: 100
均值: 49.8323
中位数: 49.7734
众数: 49.2435
标准差: 9.9247
方差: 98.4995
最小值: 27.2072
最大值: 74.4235
极差: 47.2163
四分位数Q1: 42.7174
四分位数Q2: 49.7734
四分位数Q3: 57.1504
四分位距: 14.4330
偏度: -0.0423
峰度: -0.2341
95%置信区间: [47.8743, 51.7903]
检验统计量: -0.1694
p值: 0.8652
自由度: 99
结论: 在显著性水平 0.05 下,不拒绝原假设
pearson相关系数: 0.8742
p值: 3.1643e-26
显著性: 在显著性水平 0.05 下,显著
相关强度和方向: 强相关,正相关
截距: -0.0616
斜率: 1.9774
R²值: 0.7642
调整R²值: 0.7619
F统计量: 331.1444
F检验p值: 3.16428e-26
斜率p值: 3.16428e-26
回归方程: Y = -0.0616 + 1.9774X
扩展建议
-
增强功能:
- 添加多元线性回归分析
- 实现时间序列分析功能
- 增加非参数检验方法
- 支持贝叶斯统计分析
-
用户界面:
- 开发命令行交互界面
- 创建图形界面(如使用 Tkinter 或 PyQt)
- 实现 Web 界面(如使用 Flask 或 Django)
-
性能优化:
- 针对大规模数据计算进行优化
- 添加并行计算支持
- 优化内存使用
-
教学辅助:
- 添加统计概念解释和公式推导
- 提供假设检验决策规则说明
- 实现交互式可视化(如动态调整参数)