python数学建模第九章

9.2

import numpy as np
from scipy.stats import norm, chi2, chisquare
import pylab as plt

n=50; k=6
a=np.loadtxt('ti9_2.txt').flatten()
mu=a.mean(); s=a.std(ddof=1)
x1=a.min(); x2=a.max()  #求最小值和最大值
x=np.linspace(14.55, 15.55, k)  #逐步凑的分点
bin=np.hstack([x1, x, x2])    #统计频数的分点
h=plt.hist(a, bin); f=h[0]    #统计数据频数及提取频数值
p1=norm.cdf(x, mu, s)  #计算各个分点分布函数的取值
p2=np.hstack([p1[0], np.diff(p1), 1-p1[-1]])  #计算各区间概率
print('各区间的频数：', f)
print('各区间概率为：', np.round(p2,4))
ex=n*p2  #计算期望频数
kf1=chisquare(f, ex, ddof=2)   #调用库函数计算统计量的值
kf2=sum(f**2/(n*p2))-n         #直接计算统计量的值
ka=chi2.ppf(0.95, k-2)         #临界值
print(kf1); print(round(kf2,4))  #输出两种方法计算的统计量

9.3

import numpy as np
import pylab as plt
import pandas as pd
import statsmodels.api as sm
from scipy.stats import f

a = pd.read_excel('ti9_3.xlsx', header=None)
a.boxplot()  #画箱线图
b = a.values.flatten()
print('上alpha分位数为：', f.ppf(0.95, 6, 63))
x = np.tile(np.arange(1,8),(10,1)).flatten()
d = {'x' : x, 'y' : b}  #x为类别变量
model = sm.formula.ols('y~C(x)', d).fit()  
anova = sm.stats.anova_lm(model)  #进行单因素方差分析
print('总的偏差平方和：', round(sum(anova.sum_sq),4))
print(anova); plt.show()

 

9.4

import numpy as np
import statsmodels.api as sm


y=np.loadtxt('ti9_4.txt').flatten()  #把表中的逗号替换为空格
x1=np.tile(np.arange(1,4),(12,1)).T.flatten()
x2=np.tile(np.hstack([np.ones(3),2*np.ones(3),3*np.ones(3),
                      4*np.ones(3)]),(3,1)).flatten()
d={'x1':x1, 'x2':x2, 'y':y}
md=sm.formula.ols('y~C(x1)*C(x2)',d).fit()
ano=sm.stats.anova_lm(md)
print(ano); print('总偏差平方和：', sum(ano.sum_sq))

9.5

import numpy as np
import statsmodels.api as sm

y=np.loadtxt('ti9_5.txt').T.flatten()  
x1=np.tile(np.hstack([np.ones(4),2*np.ones(4),3*np.ones(4)]),
           (4,1)).flatten()
x2=np.tile(np.tile([1,1,2,2],(1,3)),(4,1)).flatten()
x3=np.tile(np.tile([1,2],(1,6)),(4,1)).flatten()
d={'x1':x1, 'x2':x2, 'x3':x3, 'y':y}
md=sm.formula.ols('y~C(x1)*C(x2)*C(x3)',d).fit()
ano=sm.stats.anova_lm(md)
print(ano); print('总偏差平方和：', sum(ano.sum_sq))