网站 留言 以邮件形式,顺德龙江网站建设,网站策划是干什么的,页面设计第6章 数据处理与拟合模型
本章主要涉及到的知识点有#xff1a;
数据与大数据Python数据预处理常见的统计分析模型随机过程与随机模拟数据可视化
本章内容涉及到基础的概率论与数理统计理论#xff0c;如果对这部分内容不熟悉#xff0c;可以参考相关概率论与数理统计的…第6章 数据处理与拟合模型
本章主要涉及到的知识点有
数据与大数据Python数据预处理常见的统计分析模型随机过程与随机模拟数据可视化
本章内容涉及到基础的概率论与数理统计理论如果对这部分内容不熟悉可以参考相关概率论与数理统计的相关书籍如果您想对此基础做入门学习可以参考https://github.com/Git-Model/__init_Modeling__
6.1 数据与大数据
6.1.1 什么是数什么是数据
信息的不同形式称之为“模态”包括
数值类数据例如结构化的excel表格和SQL文件。文本类数据例如新闻报道、微博评论、餐饮点评等文字。图像类数据以一定尺寸的黑白或彩色图像在计算机内存储。音频类数据例如音乐、电话录音等。信号类数据例如地震波的波形、电磁波信号、脑电信号等。
这些形式下的数据都可以用来数学建模。
这些形式下的数据都可以用来数学建模。文本变为数据只需要做词频统计、TF-IDF、词嵌入等操作就可以“数化”图像变为数据因为它在计算机内本身就是以RGB三个通道的大矩阵在存储每个像素点的颜色信息音频和信号类数据则更为简单它本身可以视作一个波用信号与系统的概念和方法也可以对其“数化”。
而不同模态的数据往往能够联合对同一个事物去做描述例如对于狗狗而言既有它的图片又有它的叫声。能够同时利用描述同一事物的不同模态数据去进行联合建模的模型我们称之为多模态模型。
注意视频数据为什么没有单列呢因为视频本身也是由一帧帧图像在时间轴上堆叠才形成的数据。如果加上声音就更完美了。所以本质上视频数据是图像和音频的融合。
6.1.2 数据与大数据
真正的大数据可能要以TB甚至PB、ZB来衡量。
能够过G的数据体量基本都不小了但在G以内的数据说实话都不能算太大。
小数据有小数据的方法大数据有大数据的方法。拿着线性回归的那一套去拟合大数据往往效果不会很好拿着神经网络去学小数据学到的东西往往没有意义。面对不同的数据能够使用最合适的方法最为重要。而判断此方法是否合适一个根本的衡量就是数据的体量。
数据科学应该包含如图所示的研究方面 6.2 数据的预处理
6.2.1 为什么需要数据预处理 最上面一行我们称之为表头表头的每一格我们称之为字段。每个字段描述了这一列数据表示的意义。而这个表格的体量则是它有多少行多少列如果列数超过了行数的1/2就可以说是有些稀疏了如果列数是行数的3倍那它就是严重稀疏的数据。这个数据的每一列称其为一个属性有多少个属性又可以称之有多少维。
注意稀疏还有一种定义就是表格里面很多都是空白的或者全是0。
属性有离散属性和连续属性之分。例如在图中“品牌类型”这个属性只有{1,2,3}三个有限的取值我们称之为离散。当然这个取值也可以不是数字比如{汽车火车飞机}等。但例如属性“a1”它的取值是一个连续的浮点数这种属性我们称之为连续属性。连续属性和离散属性的处理方法是截然不同的。“目标客户编号”提供了每一行数据的标识信息可以根据这个编号对数据进行检索、排序等操作我们称之为“主码”也可以理解为id。
基础语法与一些简单的案例
#1Python创建一个数据框DataFrame
import pandas as pd
import numpy as np
data {animal: [cat, cat, snake, dog, dog, cat,
snake, cat, dog, dog],age: [2.5, 3, 0.5, np.nan, 5, 2, 4.5, np.nan, 7, 3],visits: [1, 3, 2, 3, 2, 3, 1, 1, 2, 1],priority: [yes, yes, no, yes, no, no, no,yes, no, no]}labels [a, b, c, d, e, f, g, h, i, j]
df pd.DataFrame(data)
#2显示该 DataFrame 及其数据相关的基本信息
df.describe()
#3返回DataFrame df 的前5列数据
df.head(5)
#4从 DataFrame df 选择标签列为 animal 和 age 的列
df[[animal, age]]
#5在 [3, 4, 8] 行中列为 [animal, age] 的数据
df.loc[[3, 4, 8], [animal, age]]
#6选择列为visits中等于3的行
df.loc[df[visits]3, :]
#7选择 age 为缺失值的行
df.loc[df[age].isna(), :]
#8选择 animal 是cat且age 小于 3 的行
df.loc[(df[animal] cat) (df[age] 3), :]
#9选择 age 在 2 到 4 之间的数据包含边界值
df.loc[(df[age]2)(df[age]4), :]
#10将 f 行的 age 改为 1.5
df.index labels
df.loc[[f], [age]] 1.5
#11对 visits 列的数据求和
df[visits].sum()
#12计算每种 animal age 的平均值
df.groupby([animal])[age].mean()
andas如何处理这些在数学建模中常见的任务 1创建pandas dataframe
df pd.DataFrame({From_To: [LoNDon_paris, MAdrid_miLAN, londON_StockhOlm,Budapest_PaRis, Brussels_londOn],FlightNumber: [10045, np.nan, 10065, np.nan, 10085],RecentDelays: [[23, 47], [], [24, 43, 87], [13], [67, 32]],Airline: [KLM(!), Air France (12), (British Airways. ),12. Air France, Swiss Air]})
df2FlightNumber列中有某些缺失值缺失值常用nan表示请在该列中添加10055与10075填充该缺失值。
df[FlightNumber] df[FlightNumber].interpolate().astype(int)3由于列From_To 代表从地点A到地点B因此可以将这列拆分成两列并赋予为列From与To。
temp df[From_To].str.split(_, expandTrue)
temp.columns [From, To]4将列From和To转化成只有首字母大写的形式。
temp[From] temp[From].str.capitalize()
temp[To] temp[To].str.capitalize()5将列From_To从df中去除并把列From和To添加到df中
df.drop(From_To, axis1, inplaceTrue)
df[[From, To]] temp6清除列中的特殊字符只留下航空公司的名字。
df[Airline] df[Airline].str.extract(r([a-zA-Z\s]), expandFalse).str.strip()7在 RecentDelays 列中值已作为列表输入到 DataFrame 中。我们希望每个第一个值在它自己的列中每个第二个值在它自己的列中依此类推。如果没有第 N 个值则该值应为 NaN。将 Series 列表展开为名为 的 DataFrame delays重命名列delay_1delay_2等等并将不需要的 RecentDelays 列替换df为delays。
delays df[RecentDelays].apply(pd.Series)
delays.columns [delay_%s % i for i in range(1, len(delays.columns)1)]
df df.drop(RecentDelays, axis1).join(delays, howleft)8将delay_i列的控制nan都填为自身的平均值。
for i in range(1, 4):df[fdelay_{i}] df[fdelay_{i}].fillna(np.mean(df[fdelay_{i}]))9在df中增加一行值与FlightNumber10085的行保持一致。
df df._append(df.loc[df[FlightNumber] 10085, :], ignore_indexTrue)10对df进行去重由于df添加了一行的值与FlightNumber10085的行一样的行因此去重时需要去掉。
df df.drop_duplicates()
6.2.3 数据的规约
数据为何需要规约因为实际应用中数据的分布可能是有偏的量纲影响和数值差异可能会比较大。规约是为了形成对数据的更高效表示学习到更好的模型。它会保留数据的原始特征但对极端值、异常值等会比较敏感。这里我们就介绍两个比较典型的规约方式min-max规约和Z-score规约。
min-max规约的表达式形如
这一操作的目的是为了消除量纲影响所有的属性都被规约到[0,1]的范围内数据的偏差不会那么大。但是如果出现异常值比如非常大的数值那么这个数据的分布是有偏的。为了对数据的分布进行规约还会使用到另一个常用的方法就是Z-score规约 本质上一列数据减去其均值再除以标准差如果这一列数据近似服从正态分布这个过程就是化为标准正态分布的过程。Z-score规约和min-max规约往往不是二者取其一有时候两个可以组合起来用。除此以外还有很多的规约方法但是这两种最为常用。若使用Python实现起来也非常简单。 6.3 常见的统计分析模型
说到大数据的分析思路学过概率论与数理统计这门课的同学都会脱口而出——统计分析。统计分析的很多内容与机器学习的相关内容大差不差但是侧重点不一样。统计分析通常在做“为什么”的任务机器学习通常在做“预测未来”的任务。这里将介绍几种常用的统计分析的方法由于统计分析的原理涉及大量数学知识对这部分知识的原理感兴趣的话欢迎查看https://github.com/Git-Model/Modeling-Universe/tree/main/Data-Story
6.3.1 回归分析与分类分析
回归分析中因变量是连续变量如工资、销售额而分类分析中因变量是属性变量如判断邮件“是or否”为垃圾邮件。
1回归分析在学生刚入大学时大学通常会举办一个大学入学考试该成绩虽然没什么用但是能让刚入学的学生保持持续学习的警示。现在我们想分析大学成绩GPA是否与入学成绩以及高考成绩有关
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import stats
import seaborn as sns
from IPython.display import display
import statsmodels.api as sm# 加载数据
gpa1pd.read_stata(./data/gpa1.dta)# 在数据集中提取自变量ACT为入学考试成绩、hsGPA为高考成绩
X1gpa1.ACT
X2gpa1[[ACT,hsGPA]]
# 提取因变量
ygpa1.colGPA# 为自变量增添截距项
X1sm.add_constant(X1)
X2sm.add_constant(X2)
display(X2)# 拟合两个模型
gpa_lm1sm.OLS(y,X1).fit()
gpa_lm2sm.OLS(y,X2).fit()# 输出两个模型的系数与对应p值
p1pd.DataFrame(gpa_lm1.pvalues,columns[pvalue])
c1pd.DataFrame(gpa_lm1.params,columns[params])
p2pd.DataFrame(gpa_lm2.pvalues,columns[pvalue])
c2pd.DataFrame(gpa_lm2.params,columns[params])
display(c1.join(p1,howright))
display(c2.join(p2,howright)) 从模型I可以看到ACT变量的pvalue为0.010.05说明大学入学考试成绩ACT能显著影响大学GPA但从模型II可以看到高考成绩haGPA的pvalue为0.000005远小于0.05说明高考成绩haGPA显著影响大学成绩GPA但大学入学考试成绩ACT的pvalue为0.38大于0.05不能说明大学入学考试成绩ACT显著影响大学成绩GPA。为什么会出现这样矛盾的结论呢原因有可能是在没有加入高考成绩haGPA时大学入学考试成绩ACT确实能显著影响大学成绩GPA但是加入高考成绩haGPA后有可能高考成绩haGPA高的同学大学入学考试也高分这两者的相关性较高显得“大学入学考试成绩ACT并不是那么重要”了。
从模型II的系数params数值可以看到当高考成绩haGPA显著影响大学成绩GPAparams为0.45代表每1单位的高考成绩的增加将导致大学GPA平均增加0.45个单位。
2分类分析ST是我国股市特有的一项保护投资者利于的决策当上市公司因财务状况不佳导致投资者难以预测其前景时交易所会标记该公司股票为ST并进行一系列限制措施。我们想研究被ST的公司其背后的因素并尝试通过利用公司的财务指标提前预测某上市公司在未来是否会被ST是一件很有意义的举措。而在这项任务中因变量就是公司是否会被ST。该例中自变量是一些财务指标如ARA、ASSET等它们具体的意义在此不做赘述。
# 加载基础包
import pandas as pd
import numpy as np
import statsmodels.api as sm
from scipy import stats# 读取数据
STpd.read_csv(./data/ST.csv)
ST.head()st_logitsm.formula.logit(ST~ARAASSETATOROAGROWTHLEVSHARE,
dataST).fit()
print(st_logit.summary()) 其中pvalue为第四列标记为P(| z |)的数值。与回归分析累死一般认为p值小于0.05的自变量是显著的也就是说有统计学证据表明该变量会影响因变量为1的概率即顾客购买杂志。
6.3.2 假设检验
根据样本信息与已知信息对一个描述总体性质的命题进行“是或否”的检验与回答就是假设检验的本质。即假设检验验证的不是样本本身的性质而是样本所在总体的性质。
假设检验大体上可分为两种参数假设检验与非参数假设检验。若假设是关于总体的一个参数或者参数的集合则该假设检验就是参数假设检验刚刚的例子的假设是关于总体均值的均值是参数因此这是一个参数假设检验若假设不能用一个参数集合来表示则该假设检验就是非参数假设检验一个典型就是正态性检验。
1正态性检验于参数检验比非参数检验更灵敏因此一旦数据是正态分布的我们应该使用参数检验此时对数据进行正态性检验就非常有必要了。在这里我们提供了三种方法对帮助大家判断数据的正态性可视化判断-正态分布概率图Shapiro-Wilk检验DAgostinos K-squared检验。
# 生成1000个服从正态分布的数据
data_norm stats.norm.rvs(loc10, scale10, size1000) # rvs(loc,scale,size):生成服从指定分布的随机数loc期望scale标准差size数据个数
# 生成1000个服从卡方分布的数据
data_chistats.chi2.rvs(2,3,size1000)# 画出两个概率图
figplt.figure(figsize(12,6))
ax1fig.add_subplot(1,2,1)
plot1stats.probplot(data_norm,plotax1) # 正态数据
ax2fig.add_subplot(1,2,2)
plot2stats.probplot(data_chi,plotax2) # 卡方分布数据 可以看到正态性数据在正态分布概率图中十分接近直线yx而卡方分布数据则几乎完全偏离了该直线。可见概率图可视化确实可以起到初步判断数据正态性的作用。实际应用中由于数据的复杂性仅使用一种方法判断正态性有可能产生一定的误差因此我们通常同时使用多种方法进行判断。如果不同方法得出的结论不同此时就需要仔细观察数据的特征寻找结果不一致的原因。如若Shapiro-Wilk test显著非正态DAgostinos K-squared test不显著正态则有可能是因为样本量较大或者样本中存在重复值现象如果事实确实如此那么我们就应该采纳DAgostinos K-squared test的结论而非Shapiro-Wilk test的结论。在这里我们假设数据服从正态分布p0.05则拒绝假设换个说法p0.05数据不服从正态分布如果p0.05说明没有证据说明数据不服从正态分布
# 接下来我们在python中定义一个函数将概率图、Shapiro-Wilk test、DAgostinos K-squared test结合在一起。
data_small stats.norm.rvs(0, 1, size30) # 小样本正态性数据集
data_large stats.norm.rvs(0, 1, size6000) # 大样本正态性数据集
# 定义一个正态性检验函数它可以输出
## 正态概率图
## 小样本Shapiro-Wilk检验的p值
## 大样本DAgostinos K-squared检验的p值from statsmodels.stats.diagnostic import lilliefors
from typing import Listdef check_normality(data: np.ndarray, show_flag: boolTrue) - List[float]:输入参数----------data : numpy数组或者pandas.Seriesshow_flag : 是否显示概率图Returns-------两种检验的p值概率图if show_flag:_ stats.probplot(data, plotplt)plt.show()pVals pd.Series(dtypefloat64)# DAgostinos K-squared test_, pVals[Omnibus] stats.normaltest(data) # Shapiro-Wilk test_, pVals[Shapiro-Wilk] stats.shapiro(data)print(f数据量为{len(data)}的数据集正态性假设检验的结果 : ----------------)print(pVals)
check_normality(data_small,show_flagTrue) 数据量为30的数据集正态性假设检验的结果 : ---------------- Omnibus 0.703004 Shapiro-Wilk 0.898077 dtype: float64
check_normality(data_large,show_flagFalse) # 当样本量大于5000会出现警告 数据量为6000的数据集正态性假设检验的结果 : ---------------- Omnibus 0.980014 Shapiro-Wilk 0.992427 dtype: float64
e:\anaconda\envs\ml\lib\site-packages\scipy\stats\morestats.py:1760: UserWarning: p-value may not be accurate for N 5000. warnings.warn(p-value may not be accurate for N 5000.)
2单组样本均值假定的检验必胜中学里陈老师班结束了一次英语考试。由于班级人数较多短时间内很难完成批改与统计陈老师又很想知道此次班级平均分与级长定的班级均分137的目标是否有显著区别于是他随机抽取了已经改好的10名同学的英语成绩
136,136,134,136,131,133,142,145,137,140
问陈老师可以认为此次班级平均分与级长定的班级均分137的目标没有显著区别吗
很明显这就是一个典型的单组样本均值假定的检验比较的是这个样本10个同学的英语成绩所代表的总体均值班级英语成绩均值是否与参考值137相等。那么对于这类问题我们有两种检验可以使用单样本t检验与wilcoxon检验。
定义一个单组样本均值检验函数使它可以同时输出t检验与wilcoxon符号秩和检验的p值 def check_mean(data,checkvalue,significance0.05,alternativetwo-sided): 输入参数----------data : numpy数组或者pandas.Seriescheckvalue : 想要比较的均值significance : 显著性水平alternative : 检验类型这取决于我们备择假设的符号:two-sided为双侧检验、greater为右侧检验、less为左侧检验输出-------在两种检验下的p值在显著性水平下是否拒绝原假设pValpd.Series(dtypefloat64)# 正态性数据检验-t检验_, pVal[t-test] stats.ttest_1samp(data, checkvalue,alternativealternative)print(t-test------------------------)if pVal[t-test] significance:print((目标值{0:4.2f}在显著性水平{1:}下不等于样本均值(p{2:5.3f})..format(checkvalue,significance,pVal[t-test])))else:print((目标值{0:4.2f}在显著性水平{1:}下无法拒绝等于样本均值的假设.(p{2:5.3f}).format(checkvalue,significance,pVal[t-test])))# 非正态性数据检验-wilcoxon检验_, pVal[wilcoxon] stats.wilcoxon(data-checkvalue,alternativealternative)print(wilcoxon------------------------) if pVal[wilcoxon] significance:print((目标值{0:4.2f}在显著性水平{1:}下不等于样本均值(p{2:5.3f})..format(checkvalue,significance,pVal[wilcoxon])))else:print((目标值{0:4.2f}在显著性水平{1:}下无法拒绝等于样本均值的假设.(p{2:5.3f}).format(checkvalue,significance,pVal[wilcoxon])))return pVal
2两组样本的均值相等性检验
在进行两组之间的均值比较之前我们需要进行一个重要的判断这两组样本之间是否独立呢这个问题的答案将决定着我们使用何种类型的检验。
为了让大家对这里的“独立”有较好的了解我们看看下面的例子。陈老师在年级有一个竞争对手王老师。王老师在他们班级在同一时间也举行了一次英语考试且两个班用的是同一份卷子因此两个班的英语成绩就具有比较的意义了。和陈老师一样王老师也来不及批改和统计他们班级的英语成绩不过他手头上也有12份已经改好的试卷成绩分别为
134,136,135,145,147,140,142,137,139,140,141,135
问我们可以认为两个班级的均分是相等的吗
这是一个非常典型的双独立样本的均值检验。值得注意的是这里的独立指的是抽样意义上的独立而不是统计意义的独立什么意思呢即我们只需要保证这两个样本在选取的时候是“现实上”的互不影响就可以了。至于两者在数值上是否独立通过独立性检验判断的独立性我们并不关心。对于抽样意义上的独立的解释各种教材众说纷纭在这里我们选取一种容易理解的说法两个样本中一个样本中的受试不能影响另一个样本中的受试。
在本例中两个班级的授课老师不同因此两个班学生的英语学习不会受到同一个老师的影响两个班级考试同时进行意味着不存在时间差让先考完的同学给后考完的同学通风报信进而影响受试。因此我们可以认为这是两个独立的样本。
事实上两个样本是否在抽样意义上独立是没有固定答案的因为在很多情况下我们既不能保证两个样本间完全独立也很难判断出两者是否存在相关性。我们只能在抽样的时候使用更加科学的抽样方法尽可能地避免样本间的相互影响。
若两个样本的总体都服从正态分布那么我们可以使用双样本t检验。如果不服从正态分布则可以使用Mannwhitneyu秩和检验Mannwhitneyu秩和检验是一种非参数检验。
# 定义一个单组样本均值检验函数使它可以同时输出t检验与mannwhitneyu检验的p值
def unpaired_data(group1:np.ndarray,group2:np.ndarray,significance0.05,alternativetwo-sided): 输入参数----------group1/2 : 用于比较的两组数据significance : 显著性水平alternative : 检验类型这取决于我们备择假设的符号:two-sided为双侧检验、greater为右侧检验、less为左侧检验输出-------在两种检验下的p值在显著性水平下是否拒绝原假设pValpd.Series(dtypefloat64)# 先进行两组数据的方差齐性检验_,pVal[levene]stats.levene(group1,group2)# t检验-若数据服从正态分布if pVal[levene]significance:print(在显著性水平{0:}下两组样本的方差不相等(p{1:.4f})因此需要使用方差不等的t检验.format(significance,pVal[levene]))print(------------------------------------)_, pVal[t-test] stats.ttest_ind(group1, group2,equal_varFalse,alternativealternative) # 因为方差不相等因此是Falseprint(t检验p值{}.format(pVal[t]))else:print(在显著性水平{0:}下不能拒绝两组样本方差相等的假设(p{1:.4f})因此需要使用方差相等的t检验.format(significance,pVal[levene]))print(------------------------------------)_, pVal[t-test] stats.ttest_ind(group1, group2,equal_varTrue,alternativealternative) # 因为方差相等因此是True print(t检验p值{:.3f}.format(pVal[t-test])) # mannwhitneyu检验-数据不服从正态检验_, pVal[mannwhitneyu] stats.mannwhitneyu(group1, group2,alternativealternative)print(Mann-Whitney检验p值{:.3f}.format(pVal[mannwhitneyu]))# --- STOP stats ---# 两组样本均值的散点图可视化print(------------------------------------)print(两组样本均值的散点图可视化)plt.plot(group1, bx, labelgroup1)plt.plot(group2, ro, labelgroup2)plt.legend(loc0)plt.show()return pVal
# 陈老师班
group1data
# 王老师班
group2np.array([134,136,135,145,147,140,142,137,139,140,141,135])unpaired_data(group1,group2)
在显著性水平0.05下不能拒绝两组样本方差相等的假设(p0.8277)因此需要使用方差相等的t检验
------------------------------------
t检验p值0.221
Mann-Whitney检验p值0.260
------------------------------------
两组样本均值的散点图可视化 levene 0.827727
t-test 0.221138
mannwhitneyu 0.259739
dtype: float64 两个检验都表明两个班级的平均分没有显著差异。
在进行两组间均值比较的时候有一种特殊情况——两个样本“故意”不独立。这种情况多出现两个样本分别为同一个受试个体不同时间的受试结果这两个样本是“成对”的是彼此紧密相连的。对这样两个样本进行均值比较检验就是成对检验。
为了让大家更好地理解“成对”的概念我们再把王老师的故事讲长一点。得知了王老师班的考试均值与自己班的没有显著差异陈老师非常生气“我堂堂优秀教师居然不能和隔壁老王拉开差距岂有此理”于是陈老师开始了为期一周的魔鬼训练并在一个星期后又进行了一次全班的测验。陈老师依旧抽取了上次十位同学的成绩依次为
139,141,137,136,135,132,141,148,145,139
问这次班级均分与上次是否存在显著差异呢
显然两个样本分别为相同的同学前后两次的考试成绩是非常典型的成对数据因此我们可以使用成对检验。成对检验也分为两种若总体服从正态分布则使用成对t检验若总体不服从正态分布则使用成对wilcoxon秩和检验。
def paired_data(group1:np.ndarray,group2:np.ndarray,significance,alternativetwo-sided):输入参数----------group1/2 : 用于比较的两组数据注意两组数据的样本顺序必须相同significance : 显著性水平alternative : 检验类型这取决于我们备择假设的符号:two-sided为双侧检验、greater为右侧检验、less为左侧检验输出-------在两种检验下的p值在显著性水平下是否拒绝原假设pValpd.Series(dtypefloat64)# 配对t检验-样本服从正态分布_, pVal[t-test] stats.ttest_1samp(post - pre, 0,alternativealternative)print(t-test------------------------)if pVal[t-test] significance:print((在显著性水平{0:}下两组配对样本的均值不相等(p{1:5.3f})..format(significance,pVal[t-test])))else:print((在显著性水平{0:}下无法拒绝等于样本均值的假设.(p{1:5.3f}).format(significance,pVal[t-test]))) # wilcoxon秩和检验_, pVal[wilcoxon] stats.wilcoxon(group1,group2, modeapprox,alternativealternative)print(wilcoxon------------------------)if pVal[wilcoxon] significance:print((在显著性水平{0:}下两组配对样本的均值不相等(p{1:5.3f})..format(significance,pVal[wilcoxon])))else:print((在显著性水平{0:}下无法拒绝等于样本均值的假设.(p{1:5.3f}).format(significance,pVal[wilcoxon]))) return pVal
# 第一次测验
predata
# 第二次测验
postnp.array([139,141,137,136,135,132,141,148,145,139])paired_data(pre,post,0.05)
t-test------------------------
在显著性水平0.05下两组配对样本的均值不相等(p0.039).
wilcoxon------------------------
在显著性水平0.05下两组配对样本的均值不相等(p0.049)./Users/leo/opt/anaconda3/envs/DataAnalysis_env/lib/python3.10/site-packages/scipy/stats/_morestats.py:3351: UserWarning: Sample size too small for normal approximation.warnings.warn(Sample size too small for normal approximation.)
t-test 0.039370
wilcoxon 0.048997
dtype: float64
两种检验都显示两次测验的均分在显著性水平0.05下有显著差异根据双边检验的p值是某侧单边检验两倍的结论我们可以推测出第二次测验的均分均值显著地高于第一次测验均分恭喜陈老师赢
3方差分析-多组样本间的均值相等性检验
样本间均值的“差异程度”是一个很好的评判指标但这并不足够还有一个不起眼的指标也十分重要各样本的样本内差异程度。在相同的样本间差异程度下样本内差异程度越大各总体间均值存在差异的可能性就越小。
我们需要综合这两个评判指标。最简单的方法就是两者相除即样本间均值的“差异程度”除以样本内差异程度。这就是方差分析最根本的思想。
尽管在大样本下非正态性数据的方差分析也是稳健的但是在小样本下对非正态性数据做方差分析还是可能存在误差。此时我们可以使用kruskalwallis检验。该检验也是一种非参数检验关于该检验的原理我们就不再学习了。大家只需要知道它的应用场景既可。 案例对altman_910.txt数据集进行方差分析。该数据记录了3组心脏搭桥病人给予不同水平的一氧化氮通气下他们的红细胞内叶酸水平。注三组样本都是正态性样本。
data np.genfromtxt(./data/altman_910.txt, delimiter,)
group1 data[data[:,1]1,0]
group2 data[data[:,1]2,0]
group3 data[data[:,1]3,0]
group1from typing import Tupledef anova_oneway() - Tuple[float, float]:pValpd.Series(dtypefloat64)# 先做方差齐性检验_,pVal[levene] stats.levene(group1, group2, group3)if pVal[levene]0.05: #这里假设显著性水平为0.05print(警告: 方差齐性检验的p值小于0.05: p{}方差分析结果在小样本下可能不准确.format(pVal[levene]))print(-------------------------------)# 单因素方差分析-假设样本服从正态分布_, pVal[anova_oneway_normal] stats.f_oneway(group1, group2, group3) # 在这里输入待分析的数据print(若样本服从正态分布单因素方差分析的p值为{}.format(pVal[anova_oneway_normal]))if pVal[anova_oneway_normal] 0.05:print(检验在0.05的显著性水平下显著多组样本中至少存在一组样本均值与其它样本的均值不相等。)print(---------------------------------)# 单因素方差分析-假设样本不服从正态分布_, pVal[anova_oneway_notnormal] stats.mstats.kruskalwallis(group1, group2, group3) # 在这里输入待分析的数据print(若样本不服从正态分布单因素方差分析的p值为{}.format(pVal[anova_oneway_notnormal]))if pVal[anova_oneway_notnormal] 0.05:print(检验在0.05的显著性水平下显著多组样本中至少存在一组样本均值与其它样本的均值不相等。) return pValanova_oneway()
警告: 方差齐性检验的p值小于0.05: p0.045846812634186246方差分析结果在小样本下可能不准确
-------------------------------
若样本服从正态分布单因素方差分析的p值为0.043589334959178244
检验在0.05的显著性水平下显著多组样本中至少存在一组样本均值与其它样本的均值不相等。
---------------------------------
若样本不服从正态分布单因素方差分析的p值为0.12336326887166982
levene 0.045847
anova_oneway_normal 0.043589
anova_oneway_notnormal 0.123363
dtype: float64
案例1为了改善道路的路面情况道路经常维修坑坑洼洼因此想统计一天中有多少车辆经过因为每天的车辆数都是随机的一般来说有两种技术解决这个问题
1 在道路附近安装一个计数器或安排一个技术人员在一段长时间的天数如365天每天24h统计通过道路的车辆数。
2 使用仿真技术大致模拟下道路口的场景得出一个近似可用的仿真统计指标。
由于方案1需要花费大量的人力物力以及需要花费大量的调研时间虽然能得出准确的结果但是有时候在工程应用中并不允许。因此我们选择方案2我们通过一周的简单调查得到每天的每个小时平均车辆数[30, 20, 10, 6, 8, 20, 40, 100, 250, 200, 100, 65, 100, 120, 100, 120, 200, 220, 240, 180, 150, 100, 50, 40]通过利用平均车辆数进行仿真。
# 模拟仿真研究该道路口一天平均有多少车经过
import simpyclass Road_Crossing:def __init__(self, env):self.road_crossing_container simpy.Container(env, capacity 1e8, init 0)def come_across(env, road_crossing, lmd):while True:body_time np.random.exponential(1.0/(lmd/60)) # 经过指数分布的时间后泊松过程记录数1yield env.timeout(body_time) # 经过body_time个时间yield road_crossing.road_crossing_container.put(1)hours 24 # 一天24h
minutes 60 # 一个小时60min
days 3 # 模拟3天
lmd_ls [30, 20, 10, 6, 8, 20, 40, 100, 250, 200, 100, 65, 100, 120, 100, 120, 200, 220, 240, 180, 150, 100, 50, 40] # 每隔小时平均通过车辆数
car_sum [] # 存储每一天的通过路口的车辆数之和
print(仿真开始)
for day in range(days):day_car_sum 0 # 记录每天的通过车辆数之和for hour, lmd in enumerate(lmd_ls):env simpy.Environment()road_crossing Road_Crossing(env)come_across_process env.process(come_across(env, road_crossing, lmd))env.run(until 60) # 每次仿真60minif hour % 4 0:print(第str(day1)天第str(hour1)时的车辆数, road_crossing.road_crossing_container.level)day_car_sum road_crossing.road_crossing_container.levelcar_sum.append(day_car_sum)
print(每天通过交通路口的的车辆数之和为, car_sum) 案例2现在我们来仿真“每天的商店营业额”这个复合泊松过程吧。首先我们假设每个小时进入商店的平均人数为[10, 5, 3, 6, 8, 10, 20, 40, 100, 80, 40, 50, 100, 120, 30, 30, 60, 80, 100, 150, 70, 20, 20, 10]每位顾客的平均花费为10元大约一份早餐吧请问每天商店的营业额是多少
# 模拟仿真研究该商店一天的营业额
import simpyclass Store_Money:def __init__(self, env):self.store_money_container simpy.Container(env, capacity 1e8, init 0)def buy(env, store_money, lmd, avg_money):while True:body_time np.random.exponential(1.0/(lmd/60)) # 经过指数分布的时间后泊松过程记录数1yield env.timeout(body_time) money np.random.poisson(lamavg_money)yield store_money.store_money_container.put(money)hours 24 # 一天24h
minutes 60 # 一个小时60min
days 3 # 模拟3天
avg_money 10
lmd_ls [10, 5, 3, 6, 8, 10, 20, 40, 100, 80, 40, 50, 100, 120, 30, 30, 60, 80, 100, 150, 70, 20, 20, 10] # 每个小时平均进入商店的人数
money_sum [] # 存储每一天的商店营业额总和
print(仿真开始)
for day in range(days):day_money_sum 0 # 记录每天的营业额之和for hour, lmd in enumerate(lmd_ls):env simpy.Environment()store_money Store_Money(env)store_money_process env.process(buy(env, store_money, lmd, avg_money))env.run(until 60) # 每次仿真60minif hour % 4 0:print(第str(day1)天第str(hour1)时的营业额, store_money.store_money_container.level)day_money_sum store_money.store_money_container.levelmoney_sum.append(day_money_sum)
print(每天商店的的营业额之和为, money_sum)
仿真开始
第1天第1时的营业额 121
第1天第5时的营业额 92
第1天第9时的营业额 926
第1天第13时的营业额 1036
第1天第17时的营业额 714
第1天第21时的营业额 754
第2天第1时的营业额 127
第2天第5时的营业额 123
第2天第9时的营业额 925
第2天第13时的营业额 1036
第2天第17时的营业额 567
第2天第21时的营业额 782
第3天第1时的营业额 148
第3天第5时的营业额 53
第3天第9时的营业额 1057
第3天第13时的营业额 1212
第3天第17时的营业额 605
第3天第21时的营业额 572
每天商店的的营业额之和为 [11394, 12285, 12223] 案例3艾滋病发展过程分为四个阶段状态急性感染期状态 1、无症状期状态 2 艾滋病前期状态 3, 典型艾滋病期状态 4。艾滋病发展过程基本上是一个不可逆的过程,即状态1 - 状态2 - 状态3 - 状态4。现在收集某地600例艾滋病防控数据得到以下表格 现在我们希望计算若一个人此时是无症状期状态2在10次转移之后这个人的各状态的概率是多少
# 研究无症状期病人在10期转移后的状态分布
def get_years_dist(p0, P, N):P1 Pfor i in range(N):P1 np.matmul(P1, P)return np.matmul(p0, P1)p0 np.array([0, 1, 0, 0])
P np.array([[10.0/80, 62.0/80, 5.0/80, 3.0/80],[0, 140.0/290, 93.0/290, 57.0/290],[0, 0, 180.0/220, 40.0/220],[0, 0, 0, 1]
])
N 10
print(str(N)期转移后状态分布为, np.round(get_years_dist(p0, P, N), 4)) 10期转移后状态分布为 [0.000e00 3.000e-04 1.048e-01 8.948e-01] 6.4 数据可视化
一图抵千言
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
%matplotlib inline
plt.rcParams[font.sans-serif][SimHei,Songti SC,STFangsong]
plt.rcParams[axes.unicode_minus] False # 用来正常显示负号
import seaborn as snsdf pd.DataFrame({variable: [gender, gender, age, age, age, income, income, income, income],category: [Female, Male, 1-24, 25-54, 55, Lo, Lo-Med, Med, High],value: [60, 40, 50, 30, 20, 10, 25, 25, 40],
})
df[variable] pd.Categorical(df[variable], categories[gender, age, income])
df[category] pd.Categorical(df[category], categoriesdf[category])
# 堆叠柱状图
from plotnine import *
(ggplot(df, aes(xvariable, yvalue, fillcategory))geom_col()
) # 柱状图
from plotnine import *
(ggplot(df, aes(xvariable, yvalue, fillcategory))geom_col(statidentity, positiondodge)
) 大家现在可以先忽略可视化的代码部分着重观察两个图表在表达信息上的差异。以上两个图其实都是表示在gender、age和income上不同类别的占比但是第一个图明显比第二个图表达的信息更多因为第一个图能表达每个类别的大致占比情况而第二个图只能比较不同类别下的数量情况。
6.4.1 Python三大数据可视化工具库的简介
1Matplotlib
Matplotlib正如其名脱胎于著名的建模软件Matlab因此它的设计与Matlab非常相似提供了一整套和Matlab相似的命令API适合交互式制图还可以将它作为绘图控件嵌入其它应用程序中。同时Matplotlib是Python数据可视化工具库的开山鼻祖。
Matplotlib是一个面向对象的绘图工具库pyplot是Matplotlib最常用的一个绘图接口调用方法如下
import matplotlib.pyplot as plt
在Matplotlib中我们可以想像自己手里拿着一支画笔️每一句代码都是往纸上添加一个绘图特征下面我们以最简单的方式绘制散点图为例
创建一个图形对象并设置图形对象的大小可以想象成在白纸中添加一个图并设置图的大小plt.figure(figsize(6,4))在纸上的坐标系中绘制散点plt.scatter(xx, yy)设置x轴的标签labelplt.xlabel(x)设置y轴标签的labelplt.ylabel(y)设置图表的标题plt.title(y sin(x))展示图标plt.show() 举个例子
# 创建数据
x np.linspace(-2*np.pi, 2*np.pi, 100)
y np.sin(x)
import matplotlib.pyplot as plt
plt.figure(figsize(8,6))
plt.scatter(xx, yy)
plt.xlabel(x)
plt.ylabel(y)
plt.title(y sin(x))
plt.show() 在上面的例子中我们听过Matplotlib绘制了最简单的散点图但是以上的方法没有体现Matplotlib的“面向对象”的特性。下面我们使用一个例子体会Matplotlib的面向对象绘图的特性
【例子】绘制y sin(x) 和 ycos(x)的散点图
创建第一个绘图对象 在纸上的坐标系中绘制散点plt.scatter(xx, ysin(x))设置x轴的标签labelplt.xlabel(x)设置y轴标签的labelplt.ylabel(y)设置图表的标题plt.title(y sin(x)) 创建第二个绘图对象 在纸上的坐标系中绘制散点plt.scatter(xx, ycos(x))设置x轴的标签labelplt.xlabel(x)设置y轴标签的labelplt.ylabel(y)设置图表的标题plt.title(y cos(x))
# 准备数据
x np.linspace(-2*np.pi, 2*np.pi, 100)
y1 np.sin(x)
y2 np.cos(x)
# 绘制第一个图
fig1 plt.figure(figsize(6,4), numfirst)
fig1.suptitle(y sin(x))
plt.scatter(xx, yy1)
plt.xlabel(x)
plt.ylabel(y)# 绘制第二个图
fig2 plt.figure(figsize(6,4), numsecond)
fig2.suptitle(y cos(x))
plt.scatter(xx, yy2)
plt.xlabel(x)
plt.ylabel(y)plt.show() 现在大家应该能体会到“Matplotlib的每一句代码都是往纸上添加一个绘图特征”这句话的含义了吧。现在我们来看看这样的Matplotlib有什么优点与缺点。优点是非常简单易懂而且能绘制复杂图表缺点也是十分明显的如果绘制复杂图表的时候一步一步地绘制代码量还是十分巨大的。Seaborn是在Matplotlib的基础上的再次封装是对Matplotlib绘制统计图表的简化。下面我们一起看看Seaborn的基本绘图逻辑。
2Seaborn
Seaborn主要用于统计分析绘图的它是基于Matplotlib进行了更高级的API封装。Seaborn比matplotlib更加易用尤其在统计图表的绘制上因为它避免了matplotlib中多种参数的设置。Seaborn与matplotlib关系可以把Seaborn视为matplotlib的补充。 下面我们使用一个简单的例子来看看使用Seaborn绘图与使用Matplotlib绘图之间的代码有什么不一样
# 准备数据
x np.linspace(-10, 10, 100)
y 2 * x 1 np.random.randn(100)
df pd.DataFrame({x:x, y:y})
# 准备数据
x np.linspace(-10, 10, 100)
y 2 * x 1 np.random.randn(100)
df pd.DataFrame({x:x, y:y})# 使用Seaborn绘制带有拟合直线效果的散点图
sns.lmplot(x,y,datadf) 合等统计属性高度集成在绘图函数中绘图的功能还是构筑在Matplotlib之上。因此Seaborn相当于是完善了统计图表的Matplotlib工具库二者应该是相辅相成的。因此在实际的可视化中我们往往一起使用Matplotlib和Seaborn两者的结合应该属于Python的数据可视化的一大流派吧。
3Plotnine
gplot2奠定了R语言数据可视化在R语言数据科学的统治地位R语言的数据可视化是大一统的提到R语言数据可视化首先想到的就是ggplot2。数据可视化一直是Python的短板即使有Matplotlib、Seaborn等数据可视化包也无法与R语言的ggplot2相媲美原因在于当绘制复杂图表时Matplotlib和Seaborn由于“每一句代码都是往纸上添加一个绘图特征”的特性而需要大量代码语句。Plotnine可以说是ggplot2在Python上的移植版使用的基本语法与R语言ggplot2语法“一模一样”使得Python的数据可视化能力大幅度提升为什么ggplot2和Plotnine会更适合数据可视化呢原因可以类似于PhotoShop绘图和PPT绘图的区别与PPT一笔一画的绘图方式不同的是PhotoShop绘图采用了“图层”的概念每一句代码都是相当于往图像中添加一个图层一个图层就是一类绘图动作这无疑给数据可视化工作大大减负同时更符合绘图者的认知。
下面我们通过一个案例来看看Plotnine的图层概念以及Plotnine的基本绘图逻辑
from plotnine import * # 讲Plotnine所有模块引入
from plotnine.data import * # 引入PLotnine自带数据集
mpg.head()
mpg数据集记录了美国1999年和2008年部分汽车的制造厂商型号类别驱动程序和耗油量。
manufacturer 生产厂家
model 型号类型
year 生产年份
cty 和 hwy分别记录城市和高速公路驾驶耗油量
cyl 气缸数
displ 表示发动机排量
drv 表示驱动系统前轮驱动、(f),后轮驱动®和四轮驱动(4)
class 表示车辆类型如双座汽车suv小型汽车
fl (fuel type)燃料类型 # 绘制汽车在不同驱动系统下发动机排量与耗油量的关系
p1 (ggplot(mpg, aes(xdispl, yhwy, colordrv)) # 设置数据映射图层数据集使用mpgx数据使用mpg[displ]y数据使用mpg[hwy]颜色映射使用mog[drv] geom_point() # 绘制散点图图层 geom_smooth(methodlm) # 绘制平滑线图层 labs(xdisplacement, yhorsepower) # 绘制x、y标签图层
)
print(p1) # 展示p1图像 从上面的案例我们可以看到Plotnine的绘图逻辑是一句话一个图层。因此在Plotnine中少量的代码就能画图非常漂亮的图表而且可以画出很多很复杂的图表就类似于PhotoShp能轻松画出十分复杂的图但是PPT需要大量时间也不一定能达到同样的效果。 那什么时候选择Matplotlib、Seaborn还是PlotninePlotnine具有ggplot2的图层特性但是由于开发时间较晚目前这个工具包还有一些缺陷其中最大的缺陷就是没有实现除了直角坐标以外的坐标体系如极坐标。因此Plotnine无法绘制类似于饼图、环图等图表。为了解决这个问题在绘制直角坐标系的图表时我们可以使用Plotnine进行绘制当涉及极坐标图表时我们使用Matplotlib和Seaborn进行绘制。有趣的是Matplotlib具有ggplot风格可以通过设置ggplot风格绘制具有ggplot风格的图表。
plt.style.use(ggplot) #风格使用ggplot
但是值得注意的是这里所说的绘制ggplot风格是看起来像ggplot表格但是实际上Matplotlib还是不具备图层风格。
6.4.2 基本图标Quick Start
1类别型图表类别型图表一般表现为X类别下Y数值之间的比较因此类别型图表往往包括X为类别型数据、Y为数值型数据。类别型图表常常有柱状图、横向柱状图条形图、堆叠柱状图、极坐标的柱状图、词云、雷达图、桑基图等等。
## Matplotlib绘制单系列柱状图不同城市的房价对比
data pd.DataFrame({city:[深圳, 上海, 北京, 广州, 成都], house_price(w):[3.5, 4.0, 4.2, 2.1, 1.5]})fig plt.figure(figsize(10,6))
ax1 fig.add_axes([0.15,0.15,0.7,0.7]) # [left, bottom, width, height], 它表示添加到画布中的矩形区域的左下角坐标(x, y)以及宽度和高度
plt.bar(data[city], data[house_price(w)], width0.6, aligncenter, orientationvertical, label城市)x 表示x坐标数据类型为int或float类型也可以为str
height 表示柱状图的高度也就是y坐标值数据类型为int或float类型
width 表示柱状图的宽度取值在0~1之间默认为0.8
bottom 柱状图的起始位置也就是y轴的起始坐标
align 柱状图的中心位置center,lege边缘
color 柱状图颜色
edgecolor 边框颜色
linewidth 边框宽度
tick_label 下标标签
log 柱状图y周使用科学计算方法bool类型
orientation 柱状图是竖直还是水平竖直vertical水平条horizontalplt.title(不同城市的房价对比图) # 在axes1设置标题
plt.xlabel(城市) # 在axes1中设置x标签
plt.ylabel(房价/w) # 在axes1中设置y标签
plt.grid(bTrue, whichboth) # 在axes1中设置设置网格线
for i in range(len(data)):plt.text(i-0.05, data.iloc[i,][house_price(w)]0.01, data.iloc[i,][house_price(w)],fontsize13) # 添加数据注释
plt.legend()
plt.show() ## Plotnine绘制单系列柱状图不同城市的房价对比
data pd.DataFrame({city:[深圳, 上海, 北京, 广州, 成都], house_price(w):[3.5, 4.0, 4.2, 2.1, 1.5]})p_single_bar (ggplot(data, aes(xcity, yhouse_price(w), fillcity, labelhouse_price(w)))geom_bar(statidentity)labs(x城市, y房价(w), title不同城市的房价对比图)geom_text(nudge_y0.08)theme(text element_text(family Songti SC))
)
print(p_single_bar) ## Matplotlib绘制多系列柱状图不同城市在不同年份的房价对比
data pd.DataFrame({城市:[深圳, 上海, 北京, 广州, 成都, 深圳, 上海, 北京, 广州, 成都],年份:[2021,2021,2021,2021,2021,2022,2022,2022,2022,2022],房价(w):[3.5, 4.0, 4.2, 2.1, 1.5, 4.0, 4.2, 4.3, 1.6, 1.9]
})fig plt.figure(figsize(10,6))
ax1 fig.add_axes([0.15,0.15,0.7,0.7]) # [left, bottom, width, height], 它表示添加到画布中的矩形区域的左下角坐标(x, y)以及宽度和高度
plt.bar(np.arange(len(np.unique(data[城市])))-0.15, data.loc[data[年份]2021,房价(w)], width0.3, aligncenter, orientationvertical, label年份2021)
plt.bar(np.arange(len(np.unique(data[城市])))0.15, data.loc[data[年份]2022,房价(w)], width0.3, aligncenter, orientationvertical, label年份2022)
plt.title(不同城市的房价对比图) # 在axes1设置标题
plt.xlabel(城市) # 在axes1中设置x标签
plt.ylabel(房价/w) # 在axes1中设置y标签
plt.xticks(np.arange(len(np.unique(data[城市]))), np.array([深圳, 上海, 北京, 广州, 成都]))
plt.grid(bTrue, whichboth) # 在axes1中设置设置网格线data_2021 data.loc[data[年份]2021,:]
for i in range(len(data_2021)):plt.text(i-0.15-0.05, data_2021.iloc[i,2]0.05, data_2021.iloc[i,2],fontsize13) # 添加数据注释data_2022 data.loc[data[年份]2022,:]
for i in range(len(data_2022)):plt.text(i0.15-0.05, data_2022.iloc[i,2]0.05, data_2022.iloc[i,2],fontsize13) # 添加数据注释
plt.legend()
plt.show() ## Plotnine绘制多系列柱状图不同城市在不同年份的房价对比
data pd.DataFrame({城市:[深圳, 上海, 北京, 广州, 成都, 深圳, 上海, 北京, 广州, 成都],年份:[2021,2021,2021,2021,2021,2022,2022,2022,2022,2022],房价(w):[3.5, 4.0, 4.2, 2.1, 1.5, 4.0, 4.2, 4.3, 1.6, 1.9]
})data[年份] pd.Categorical(data[年份], orderedTrue, categoriesdata[年份].unique())
p_mult_bar (ggplot(data, aes(x城市, y房价(w), fill年份))geom_bar(statidentity,width0.6, positiondodge)scale_fill_manual(values [#f6e8c3, #5ab4ac])labs(x城市, y房价(w), title不同城市的房价对比图)geom_text(aes(label房价(w)), position position_dodge2(width 0.6, preserve single))theme(text element_text(family Songti SC))
)
print(p_mult_bar) ## Matplotlib绘制堆叠柱状图不同城市在不同年份的房价对比
data pd.DataFrame({城市:[深圳, 上海, 北京, 广州, 成都, 深圳, 上海, 北京, 广州, 成都],年份:[2021,2021,2021,2021,2021,2022,2022,2022,2022,2022],房价(w):[3.5, 4.0, 4.2, 2.1, 1.5, 4.0, 4.2, 4.3, 1.6, 1.9]
})
tmpdata.set_index([城市,年份])[房价(w)].unstack()
datatmp.rename_axis(columnsNone).reset_index()
data.columns [城市,2021房价,2022房价]
print(data)plt.figure(figsize(10,6))
plt.bar(data[城市], data[2021房价], width0.6, aligncenter, orientationvertical, label年份2021)
plt.bar(data[城市], data[2022房价], width0.6, aligncenter, orientationvertical, bottomdata[2021房价],label年份2022)
plt.title(不同城市2121-2022年房价对比图) # 在axes1设置标题
plt.xlabel(城市) # 在axes1中设置x标签
plt.ylabel(房价/w) # 在axes1中设置y标签
plt.legend()
plt.show() ## Matplotlib绘制百分比柱状图不同城市在不同年份的房价对比
## Matplotlib绘制堆叠柱状图不同城市在不同年份的房价对比
data pd.DataFrame({城市:[深圳, 上海, 北京, 广州, 成都, 深圳, 上海, 北京, 广州, 成都],年份:[2021,2021,2021,2021,2021,2022,2022,2022,2022,2022],房价(w):[3.5, 4.0, 4.2, 2.1, 1.5, 4.0, 4.2, 4.3, 1.6, 1.9]
})
tmpdata.set_index([城市,年份])[房价(w)].unstack()
datatmp.rename_axis(columnsNone).reset_index()
data.columns [城市,2021房价,2022房价]
print(data)plt.figure(figsize(10,6))
plt.bar(data[城市], data[2021房价]/(data[2021房价]data[2022房价]), width0.4, aligncenter, orientationvertical, label年份2021)
plt.bar(data[城市], data[2022房价]/(data[2021房价]data[2022房价]), width0.4, aligncenter, orientationvertical, bottomdata[2021房价]/(data[2021房价]data[2022房价]),label年份2022)
plt.title(不同城市2121-2022年房价对比图) # 设置标题
plt.xlabel(城市) # 在axes1中设置x标签
plt.ylabel(房价/w) # 在axes1中设置y标签
plt.legend()
plt.show() # 使用Matplotlib绘制雷达图英雄联盟几位英雄的能力对比
data pd.DataFrame({属性: [血量, 攻击力, 攻速, 物抗, 魔抗],艾希:[3, 7, 8, 2, 2],诺手:[8, 6, 3, 6, 6]
})plt.figure(figsize(8,8))
theta np.linspace(0, 2*np.pi, len(data), endpointFalse) # 每个坐标点的位置
theta np.append(theta, theta[0]) # 让数据封闭
aixi np.append(data[艾希].values,data[艾希][0]) #让数据封闭
nuoshou np.append(data[诺手].values,data[诺手][0]) # 让数据封闭
shuxing np.append(data[属性].values,data[属性][0]) # 让数据封闭plt.polar(theta, aixi, ro-, lw2, label艾希) # 画出雷达图的点和线
plt.fill(theta, aixi, facecolorr, alpha0.5) # 填充
plt.polar(theta, nuoshou, bo-, lw2, label诺手) # 画出雷达图的点和线
plt.fill(theta, nuoshou, facecolorb, alpha0.5) # 填充
plt.thetagrids(theta/(2*np.pi)*360, shuxing) # 为每个轴添加标签
plt.ylim(0,10)
plt.legend()
plt.show() 2关系型图表关系型图表一般表现为X数值与Y数值之间的关系如是否是线性关系、是否有正向相关关系等等。一般来说关系可以分为数值型关系、层次型关系和网络型关系。
# 使用Matplotlib和四个图说明相关关系
x np.random.randn(100)*10
y1 np.random.randn(100)*10
y2 2 * x 1 np.random.randn(100)
y3 -2 * x 1 np.random.randn(100)
y4 x**2 1 np.random.randn(100)plt.figure(figsize(12, 12))plt.subplot(2,2,1) #创建两行两列的子图并绘制第一个子图
plt.scatter(x, y1, cdodgerblue, marker., s50)
plt.xlabel(x)
plt.ylabel(y1)
plt.title(y1与x不存在关联关系)plt.subplot(2,2,2) #创建两行两列的子图并绘制第二个子图
plt.scatter(x, y2, ctomato, markero, s10)
plt.xlabel(x)
plt.ylabel(y2)
plt.title(y2与x存在关联关系)plt.subplot(2,2,3) #创建两行两列的子图并绘制第三个子图
plt.scatter(x, y3, cmagenta, markero, s10)
plt.xlabel(x)
plt.ylabel(y3)
plt.title(y3与x存在关联关系)plt.subplot(2,2,4) #创建两行两列的子图并绘制第四个子图
plt.scatter(x, y4, cdeeppink, markers, s10)
plt.xlabel(x)
plt.ylabel(y3)
plt.title(y4与x存在关联关系)plt.show() # 使用Plotnine和四个图说明相关关系
x np.random.randn(100)*10
y1 np.random.randn(100)*10
y2 10 * x 1 np.random.randn(100)
y3 -10 * x 1 np.random.randn(100)
y4 x**2 1 np.random.randn(100)df pd.DataFrame({x: np.concatenate([x,x,x,x]),y: np.concatenate([y1, y2, y3, y4]),class: [y1]*100 [y2]*100 [y3]*100 [y4]*100
})p1 (ggplot(df)geom_point(aes(xx, yy, fillclass, shapeclass), colorblack, size2)scale_fill_manual(values(#00AFBB, #FC4E07, #00589F, #F68F00))theme(text element_text(family Songti SC))
)
print(p1) # 使用Matplotlib绘制具备趋势线的散点图
from sklearn.linear_model import LinearRegression #线性回归等参数回归
import statsmodels.api as smfrom sklearn.preprocessing import PolynomialFeatures # 构造多项式
x np.linspace(-10, 10, 100)
y np.square(x) np.random.randn(100)*100
x_poly2 PolynomialFeatures(degree2).fit_transform(x.reshape(-1, 1))
y_linear_pred LinearRegression().fit(x.reshape(-1, 1), y).predict(x.reshape(-1, 1))
y_poly_pred LinearRegression().fit(x_poly2, y).predict(x_poly2)
y_exp_pred LinearRegression().fit(np.exp(x).reshape(-1, 1), y).predict(np.exp(x).reshape(-1, 1))
y_loess_pred sm.nonparametric.lowess(x, y, frac2/3)[:, 1]plt.figure(figsize(8, 8))
plt.subplot(2,2,1)
plt.scatter(x, y, ctomato, markero, s10)
plt.plot(x, y_linear_pred, cdodgerblue)
plt.xlabel(x)
plt.ylabel(y)
plt.title(带线性趋势线的散点图)plt.subplot(2,2,2)
plt.scatter(x, y, ctomato, markero, s10)
plt.plot(x, y_poly_pred, cdodgerblue)
plt.xlabel(x)
plt.ylabel(y)
plt.title(带二次趋势线的散点图)plt.subplot(2,2,3)
plt.scatter(x, y, ctomato, markero, s10)
plt.plot(x, y_exp_pred, cdodgerblue)
plt.xlabel(x)
plt.ylabel(y)
plt.title(带指数趋势线的散点图)plt.subplot(2,2,4)
plt.scatter(x, y, ctomato, markero, s10)
plt.plot(x, y_loess_pred, cdodgerblue)
plt.xlabel(x)
plt.ylabel(y)
plt.title(带 loess平滑线的散点图)plt.show() # 使用Matplotlib绘制聚类散点图
from sklearn.datasets import load_iris #家在鸢尾花数据集
iris load_iris()
X iris.data
label iris.target
feature iris.feature_names
df pd.DataFrame(X, columnsfeature)
df[label] labellabel_unique np.unique(df[label]).tolist()
plt.figure(figsize(10, 6))
for i in label_unique:df_label df.loc[df[label] i, :]plt.scatter(xdf_label[sepal length (cm)], ydf_label[sepal width (cm)], s20, labeli)
plt.xlabel(sepal length (cm))
plt.ylabel(sepal width (cm))
plt.title(sepal width (cm)~sepal length (cm))
plt.legend()
plt.show() # 使用plotnine绘制相关系数矩阵图
from plotnine.data import mtcars
corr_mat np.round(mtcars.corr(), 1).reset_index() #计算相关系数矩阵
df pd.melt(corr_mat, id_varsindex, var_namevariable, value_namecorr_xy) #将矩阵宽数据变成长数据
df[abs_corr] np.abs(df[corr_xy])
p1 (ggplot(df, aes(xindex, yvariable, fillcorr_xy, sizeabs_corr))geom_point(shapeo, colorblack)scale_size_area(max_size11, guideFalse)scale_fill_cmap(nameRdYIBu_r)coord_equal()labs(xVariable, yVariable)theme(dpi100, figure_size(4.5,4.55))
)
p2 (ggplot(df, aes(xindex, yvariable, fillcorr_xy, sizeabs_corr))geom_point(shapes, colorblack)scale_size_area(max_size10, guideFalse)scale_fill_cmap(nameRdYIBu_r)coord_equal()labs(xVariable, yVariable)theme(dpi100, figure_size(4.5,4.55))
)
p3 (ggplot(df, aes(xindex, yvariable, fillcorr_xy, labelcorr_xy))geom_tile(colorblack)geom_text(size8, colorwhite)scale_fill_cmap(nameRdYIBu_r)coord_equal()labs(xVariable, yVariable)theme(dpi100, figure_size(4.5,4.55))
)
print([p1, p2, p3])
# 使用plotnine绘制相关系数矩阵图
from plotnine.data import mtcars
corr_mat np.round(mtcars.corr(), 1).reset_index() #计算相关系数矩阵
df pd.melt(corr_mat, id_varsindex, var_namevariable, value_namecorr_xy) #将矩阵宽数据变成长数据
df[abs_corr] np.abs(df[corr_xy])
p1 (ggplot(df, aes(xindex, yvariable, fillcorr_xy, sizeabs_corr))geom_point(shapeo, colorblack)scale_size_area(max_size11, guideFalse)scale_fill_cmap(nameRdYIBu_r)coord_equal()labs(xVariable, yVariable)theme(dpi100, figure_size(4.5,4.55))
)
p2 (ggplot(df, aes(xindex, yvariable, fillcorr_xy, sizeabs_corr))geom_point(shapes, colorblack)scale_size_area(max_size10, guideFalse)scale_fill_cmap(nameRdYIBu_r)coord_equal()labs(xVariable, yVariable)theme(dpi100, figure_size(4.5,4.55))
)
p3 (ggplot(df, aes(xindex, yvariable, fillcorr_xy, labelcorr_xy))geom_tile(colorblack)geom_text(size8, colorwhite)scale_fill_cmap(nameRdYIBu_r)coord_equal()labs(xVariable, yVariable)theme(dpi100, figure_size(4.5,4.55))
)
print([p1, p2, p3]) # 使用Matplotlib/Seaborn绘制相关系数矩阵图
uniform_data np.random.rand(10, 12)
sns.heatmap(uniform_data) 3分布型图表所谓数据的分布其实就是数据在哪里比较密集那里比较稀疏描述数据的密集或者稀疏情况实际上可以用频率或者概率。
# 使用matplotlib绘制直方图
plt.figure(figsize(8, 6))
plt.hist(mtcars[mpg], bins20, alpha0.85)
plt.xlabel(mpg)
plt.ylabel(count)
plt.show() # 使用matplotlib绘制箱线图
import seaborn as sns
from plotnine.data import mtcarsdata mtcars
data[carb] data[carb].astype(category)
plt.figure(figsize(8, 6))
sns.boxenplot(xcarb, ympg, datamtcars, linewidth0.2, palettesns.husl_palette(6, s0.9, l0.65, h0.0417))
plt.show() # 使用Matplotlib绘制饼状图
from matplotlib import cm, colors
df pd.DataFrame({己方: [寒冰, 布隆, 发条, 盲僧, 青钢影],敌方: [女警, 拉克丝, 辛德拉, 赵信, 剑姬],己方输出: [26000, 5000, 23000, 4396, 21000],敌方输出: [25000, 12000, 21000, 10000, 18000]
})df_our df[[己方, 己方输出]].sort_values(by己方输出, ascendingFalse).reset_index()
df_other df[[敌方, 敌方输出]].sort_values(by敌方输出, ascendingFalse).reset_index()
color_list [cm.Set3(i) for i in range(len(df))]
plt.figure(figsize(16, 10))
plt.subplot(1,2,1)
plt.pie(df_our[己方输出].values, startangle90, shadowTrue, colorscolor_list, labelsdf_our[己方].tolist(), explode(0,0,0,0,0.3), autopct%.2f%%)plt.subplot(1,2,2)
plt.pie(df_other[敌方输出].values, startangle90, shadowTrue, colorscolor_list, labelsdf_other[敌方].tolist(), explode(0,0,0,0,0.3), autopct%.2f%%)plt.show() # 使用Matplotlib绘制环状图
from matplotlib import cm, colors
df pd.DataFrame({己方: [寒冰, 布隆, 发条, 盲僧, 青钢影],敌方: [女警, 拉克丝, 辛德拉, 赵信, 剑姬],己方输出: [26000, 5000, 23000, 4396, 21000],敌方输出: [25000, 12000, 21000, 10000, 18000]
})df_our df[[己方, 己方输出]].sort_values(by己方输出, ascendingFalse).reset_index()
df_other df[[敌方, 敌方输出]].sort_values(by敌方输出, ascendingFalse).reset_index()
color_list [cm.Set3(i) for i in range(len(df))]
wedgeprops {width:0.3, edgecolor:black, linewidth:3}
plt.figure(figsize(16, 10))
plt.subplot(1,2,1)
plt.pie(df_our[己方输出].values, startangle90, shadowTrue, colorscolor_list, wedgepropswedgeprops, labelsdf_our[己方].tolist(), explode(0,0,0,0,0.3), autopct%.2f%%)
plt.text(0, 0, 己方 , hacenter, vacenter, fontsize30)plt.subplot(1,2,2)
plt.pie(df_other[敌方输出].values, startangle90, shadowTrue, colorscolor_list, wedgepropswedgeprops, labelsdf_other[敌方].tolist(), explode(0,0,0,0,0.3), autopct%.2f%%)
plt.text(0, 0, 敌方 , hacenter, vacenter, fontsize30)
plt.show() 4时间序列型图表
# 使用Plotnine绘制时间序列线图
df pd.read_csv(./data/AirPassengers.csv) # 航空数据1949-1960
df[date] pd.to_datetime(df[date])
p1 (ggplot(df, aes(xdate, yvalue))geom_line(size1, colorred)scale_x_date(date_labels%Y, date_breaks1 year)xlab(date)ylab(value)
)
print(p1) # 使用Matplotlib绘制时间序列折线图
plt.figure(figsize(8,6))
plt.plot(df[date], df[value], colorred)
plt.xlabel(date)
plt.ylabel(value)
plt.show() # Plotnine绘制多系列折线图
date_list pd.date_range(2022-01-01, 2022-03-31).astype(str).tolist() * 2
value_list np.random.randn(len(date_list))
Class [1] * (len(date_list) // 2) [2] * (len(date_list) // 2)
data pd.DataFrame({date_list: date_list,value_list: value_list,Class: Class
})data[date_list] pd.to_datetime(data[date_list])
p1 (ggplot(data, aes(xdate_list, yvalue_list, groupfactor(Class), colorfactor(Class)))geom_line(size1)scale_x_date(date_labels%D)xlab(date)ylab(value)
)
print(p1) # Matplotlib 绘制多系列折线图
date_list pd.date_range(2022-01-01, 2022-03-31).astype(str).tolist()
value_list1 np.random.randn(len(date_list))
value_list2 np.random.randn(len(date_list))
data pd.DataFrame({date_list: date_list,value_list1: value_list1,value_list2: value_list2
})
data[date_list] pd.to_datetime(data[date_list])plt.figure(figsize(8, 6))
plt.plot(data[date_list], data[value_list1], colorred, alpha0.86, labelvalue1)
plt.plot(data[date_list], data[value_list2], colorblue, alpha0.86, labelvalue2)
plt.legend()
plt.xlabel(date)
plt.ylabel(value)
plt.show() 6.5 插值模型
6.5.1 线性插值法
线性插值对两个点中的解析式是按照线性方程来建模。比如我们得到的原始数据列{y}和数据的下标{x}这里数据下标x可能并不是固定频率的连续取值而是和y一样存在缺失的。给定了数据点(xk,yk)和(xk1,yk1)需要对两个点之间构造直线进行填充。很显然根据直线的点斜式方程这个直线解析式为 $$ {L_1}(x) {y_k} \frac{{{y_{k 1}} - {y_k}}}{{{x_{k 1}} - {x_k}}}(x - {x_k}) \tag{6.5.1} $$ 按照这个方程形式对空缺的(x,y)处进行填充就可以完成插值过程。 import numpy as np
#数据准备
X np.arange(-2*np.pi, 2*np.pi, 1) # 定义样本点X从-pi到pi每次间隔1
Y np.sin(X) # 定义样本点Y形成sin函数
new_x np.arange(-2*np.pi, 2*np.pi, 0.1) # 定义插值点
# 进行样条插值
import scipy.interpolate as spi
# 进行一阶样条插值
ipo1 spi.splrep(X, Y, k1) # 样本点导入生成参数
iy1 spi.splev(new_x, ipo1) # 根据观测点和样条参数生成插值
6.4.2 三次样条插值
三次样条插值是一种非常自然的插值方法。它将两个数据点之间的填充模式设置为三次多项式。它假设在数据点(xk,yk)和(xk1,yk1)之间的三次式叫做Ik那么这一组三次式需要满足条件 $$ {a_i}x_i^3 {b_i}x_i^2 {c_i}{x_i} {d_i} {a_{i 1}}x_{i 1}^3 {b_{i 1}}x_{i 1}^2 {c_{i 1}}{x_{i 1}} {d_{i 1}} \tag{6.5.2} $$ $$ 3{a_i}x_i^2 2{b_i}{x_i} {c_i} 3{a_{i 1}}x_{i 1}^2 2{b_{i 1}}{x_{i 1}} {c_{i 1}} \tag{6.5.3} $$ $$ 6{a_i}{x_i} 2{b_i} 6{a_{i 1}}{x_{i 1}} 2{b_{i 1}} \tag{6.5.4} $$ 通过解方程的形式可以解出每一只三次式。而简而言之也就是说某个数据点前后两条三次函数不仅在当前的函数值相等一次导数值和二次导数值也要保持相等。
# 进行三次样条拟合
ipo3 spi.splrep(X, Y, k3) # 样本点导入生成参数
iy3 spi.splev(new_x, ipo3) # 根据观测点和样条参数生成插值
6.4.3 拉格朗日插值
对于一组数据{y}和下标{x}定义n个拉格朗日插值基函数 $$ {l_k}(x) \prod\limits_{i 0,i \ne k}^n {\frac{{x - {x_i}}}{{{x_k} - {x_i}}}} \tag{6.5.5} $$ 这本质上是一个分式当xxk时lk(x)1这一操作实现了离散数据的连续化。按照对应下标的函数值加权求和可以得到整体的拉格朗日插值函数 $$ L(x) \sum\limits_{k 0}^n {{y_k}{l_k}(x)} \tag{6.5.6} $$ 但很不幸Python没有为我们提供拉格朗日插值函数的方法。不过好在实现并不算太困难我们可以自己写一个
def lagrange(x0,y0,x):y[]for k in range(len(x)):s0for i in range(len(y0)):ty0[i]for j in range(len(y0)):if i!j:t*(x[k]-x0[j])/(x0[i]-x0[j])sty.append(s)return y 通过一个简单的样例进行测试
from scipy.interpolate import interp1d
x0[1,2,3,4,5]
y0[1.6,1.8,3.2,5.9,6.8]
xnp.arange(1,5,1/30)
f1interp1d(x0,y0,linear)
y1f1(x)
f2interp1d(x0,y0,cubic)
y2f2(x)
y3lagrange(x0,y0,x)
plt.plot(x0,y0,r*)
plt.plot(x,y1,b-,x,y2,y-,x,y3,r-)
plt.show()
