Datawhale学习打卡-数学建模导论(下)Task1

数据处理与拟合模型

数学建模中的数据处理以及常见的数据模型

Part 1 数据与大数据

1.1 何为数据?

        数据存在是为了描述与传递多种多样的信息。其中,信息的不同形式被成为“模态”,包括:

  • 数值类数据,例如结构化的excel表格和SQL文件。
  • 文本类数据,例如新闻报道、微博评论、餐饮点评等文字。
  • 图像类数据,以一定尺寸的黑白或彩色图像在计算机内存储。
  • 音频类数据,例如音乐、电话录音等。
  • 信号类数据,例如地震波的波形、电磁波信号、脑电信号等。
1.2 如何把信息转变为数据?

        依据信息特点对信息进行处理,如:

  • 文本变为数据,只需要做词频统计、TF-IDF、词嵌入等操作就可以“数化”;
  • 图像变为数据,因为它在计算机内本身就是以RGB三个通道的大矩阵在存储每个像素点的颜色信息;
  • 音频和信号类数据本身可以视作一个波,用信号与系统的概念和方法也可以对其“数化”。
  • 视频本身是由一帧帧图像在时间轴上堆叠才形成的数据,是图像和音频的融合。

     其中,能够同时利用描述同一事物的不同模态数据去进行联合建模的模型,我们称之为多模态模型。

1.3 何为大数据?

大数据(Big Data)是指数据量巨大、类型多样、处理速度快、价值密度低的数据集合。

大数据具有以下几个显著特点:

  1. 体量大(Volume):数据量巨大,通常以TB(太字节)、PB(拍字节)甚至EB(艾字节)为单位。
  2. 多样性(Variety):数据类型多样,包括结构化数据(如数据库中的表格数据)、半结构化数据(如XML、JSON)、非结构化数据(如文本、图片、视频)等。
  3. 处理速度快(Velocity):数据的处理和分析需要快速响应,实时或近实时处理能力是大数据应用的关键。
  4. 价值密度低(Value):虽然数据量巨大,但其中有价值的信息往往只占很小的一部分,需要通过复杂的分析和挖掘技术来提取。
  5. 真实性(Veracity):数据的准确性和可靠性是大数据应用的基础,数据的质量问题会直接影响分析结果的有效性。

大数据技术的核心包括数据的获取与存储、数据处理、数据分析、管理及数据挖掘、应用等。

注:① 面对不同的数据,能够使用最合适的方法最为重要。而判断此方法是否合适,一个根本的衡量就是数据的体量。

       ② csv、tsv等格式的数据可以用记事本、excel、WPS等打开。

Part 2 Python数据预处理

2.1 数据预处理

        以表格为例,

 图 1 2021年华数杯C题的部分数据在WPS中展示图

  • 表头:最上面一行
  • 字段:表头的每一格,每个字段描述了这一列数据表示的意义
  • 体量:指表格有多少行多少列,如果列数超过了行数的1/2就可以说有些稀疏,如果列数是行数的3倍,那它就是严重稀疏的数据(稀疏还有一种定义,就是表格里面很多都是空白的或者全是0)。
  • 属性(维):数据的每一列称其为一个属性,有多少个属性又可以称之有多少维。其中,属性有离散属性和连续属性之分。例如在图中,“品牌类型”这个属性只有{1,2,3}三个有限的取值,我们称之为离散。当然,这个取值也可以不是数字,比如{汽车,火车,飞机}等。但例如属性“a1”,它的取值是一个连续的浮点数,这种属性我们称之为连续属性。连续属性和离散属性的处理方法是截然不同的。
2.2 不规则数据的处理

        对于一份非理想的不规则数据经常会出现空缺、重复和一些异常记录,需要采取一定策略对其进行识别并处理。

  • 对于重复数据直接将其删除即可。
  • 缺失数据的处理主要是观察缺失率,如果存在缺失的数据项(数表一行称作一个数据项)占比较少(大概5%以内)这个时候如果问题允许可以把行删掉。如果缺失率稍微高一点(5%-20%)左右就可以使用填充、插值的方法去处理,而填充方法包括常数填充、均值填充等方法;如果缺失率还高一些(20%-40%)那么就需要用预测方法例如机器学习去填充缺失数据了;但如果一行数据有50%以上都是缺失的,如果条件允许,我们可以把这一列都删掉(凡事都有例外,见机行事)。
2.3 使用python处理数据

2.3.1 常见的Python数据处理库及其主要功能

NumPy

  • 功能:提供多维数组对象、派生对象(如掩码数组和矩阵)以及用于快速操作数组的各种例程,包括数学、逻辑、形状操作、排序、选择、I/O、离散傅立叶变换、基本线性代数、基本统计运算、随机模拟等。
  • 用途:科学计算、数据分析、机器学习。

Pandas

  • 功能:提供高性能、易用的数据结构和数据分析工具。主要数据结构是Series(一维)和DataFrame(二维),支持时间序列分析。
  • 用途:数据清洗、处理、分析,时间序列数据。

SciPy

  • 功能:基于NumPy,提供更多的科学计算工具,如优化、线性代数、积分、插值、特殊函数、快速傅立叶变换、信号处理和统计学。
  • 用途:科学计算、数据分析。

Matplotlib

  • 功能:提供类似MATLAB的绘图框架,支持多种输出格式,能够绘制各种静态、动态、交互式的图表。
  • 用途:数据可视化。

Seaborn

  • 功能:基于Matplotlib,提供更高级的统计图表绘制功能,支持多种图表类型,如分布图、箱型图、热力图等。
  • 用途:数据可视化,统计分析。

Scikit-learn

  • 功能:提供简单有效的工具进行数据挖掘和数据分析,包括分类、回归、聚类、降维、模型选择和评估等。
  • 用途:机器学习、数据挖掘。

Statsmodels

  • 功能:提供统计模型的估计和统计测试,包括线性回归、时间序列分析、面板数据等。
  • 用途:统计分析、经济计量学。

TensorFlowPyTorch

  • 功能:提供深度学习框架,支持神经网络的构建、训练和测试。
  • 用途:深度学习、机器学习。

以一个简单的demo,说明pandas如何处理这些在数学建模中常见的任务:

#(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"']})
df

#(2)FlightNumber列中有某些缺失值,缺失值常用nan表示,请在该列中添加10055与10075填充该缺失值。
df['FlightNumber'] = df['FlightNumber'].interpolate().astype(int)

#(3)由于列From_To 代表从地点A到地点B,因此可以将这列拆分成两列,并赋予为列From与To。
temp = df['From_To'].str.split("_", expand=True)
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', axis=1, inplace=True)
df[['From', 'To']] = temp

#(6)清除列中的特殊字符,只留下航空公司的名字。
df['Airline'] = df['Airline'].str.extract(r'([a-zA-Z\s]+)', expand=False).str.strip()

#(7)在 RecentDelays 列中,值已作为列表输入到 DataFrame 中。我们希望每个第一个值在它自己的列中,每个第二个值在它自己的列中,依此类推。如果没有第 N 个值,则该值应为 NaN。将 Series 列表展开为名为 的 DataFrame delays,重命名列delay_1,delay_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', axis=1).join(delays, how='left')

#(8)将delay_i列的控制nan都填为自身的平均值。
for i in range(1, 4):
    df[f'delay_{i}'] = df[f'delay_{i}'].fillna(np.mean(df[f'delay_{i}']))

#(9)在df中增加一行,值与FlightNumber=10085的行保持一致。
df = df._append(df.loc[df['FlightNumber'] == 10085, :], ignore_index=True)

#(10)对df进行去重,由于df添加了一行的值与FlightNumber=10085的行一样的行,因此去重时需要去掉。
df = df.drop_duplicates()
2.4 数据的规约

2.4.1 数据为何需要规约?

        因为实际应用中数据的分布可能是有偏的,量纲影响和数值差异可能会比较大。规约是为了形成对数据的更高效表示,学习到更好的模型。它会保留数据的原始特征,但对极端值、异常值等会比较敏感

min-max规约

目的:消除量纲影响,所有的属性都被规约到[0,1]的范围内,数据的偏差不会那么大。但是如果出现异常值,比如非常大的数值,那么这个数据的分布是有偏的。

Z-score规约

 一列数据减去其均值再除以标准差,如果这一列数据近似服从正态分布,这个过程就是化为标准正态分布的过程。

注:Z-score规约和min-max规约往往不是二者取其一,有时候两个可以组合起来用。

Part 3 常见的统计分析模型

3.1 回归分析与分类分析

目的:研究因变量(被解释变量)与自变量(解释变量)之间存在的潜在关系,并通过统计模型的形式将这些潜在关系进行显式的表达。其中,回归分析中因变量是连续变量,如工资、销售额;分类分析中因变量是属性变量,如判断邮件“是or否”为垃圾邮件。

(1)回归分析示例——分析大学成绩GPA是否与入学成绩以及高考成绩有关

其中,gpa_lm1=sm.OLS(y,X1).fit() 是一个使用Python中的 statsmodels 库进行普通最小二乘回归分析的代码片段,指使用 statsmodels 库中的 OLS 方法,以 X1 作为自变量,y 作为因变量,拟合一个线性回归模型,并将拟合后的模型保存在变量 gpa_lm1 中。

  • sm 是 statsmodels 的缩写,是Python中用于统计和经济计量学的一个库。
  • OLS 代表“普通最小二乘法”(Ordinary Least Squares),是一种常用的线性回归方法。
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

# 加载数据
gpa1=pd.read_stata('./data/gpa1.dta')

# 在数据集中提取自变量(ACT为入学考试成绩、hsGPA为高考成绩)
X1=gpa1.ACT
X2=gpa1[['ACT','hsGPA']]
# 提取因变量
y=gpa1.colGPA

# 为自变量增添截距项
X1=sm.add_constant(X1)
X2=sm.add_constant(X2)
display(X2)


# 拟合两个模型
gpa_lm1=sm.OLS(y,X1).fit()
gpa_lm2=sm.OLS(y,X2).fit()

# 输出两个模型的系数与对应p值
p1=pd.DataFrame(gpa_lm1.pvalues,columns=['pvalue'])
c1=pd.DataFrame(gpa_lm1.params,columns=['params'])
p2=pd.DataFrame(gpa_lm2.pvalues,columns=['pvalue'])
c2=pd.DataFrame(gpa_lm2.params,columns=['params'])
display(c1.join(p1,how='right'))
display(c2.join(p2,how='right'))

样例输出结果分析:

在使用 statsmodels 进行线性回归分析时,paramspvalue 是模型拟合结果中的两个重要属性:

  1. params:是模型参数的估计值。在普通最小二乘回归中,params 包含了回归系数的估计值。对于一个简单的线性回归模型 𝑦=𝛽0+𝛽1𝑥1+𝜖,params 将包含截距项 𝛽0β0​ 和斜率 𝛽1β1​ 的估计值。例如,如果 X1 包含一个变量 x1,那么 params 将返回一个数组,其中第一个元素是截距项的估计值,第二个元素是 x1 的斜率系数的估计值。

  2. pvalue:是模型参数的 p 值。p 值用于统计检验,判断模型中的系数是否显著不同于零。在统计学中,通常使用 p 值来确定一个效应(例如,自变量对因变量的影响)是否显著。具体来说,p 值是观察到的统计量(如 t 统计量)或更极端情况下的概率。如果 p 值小于某个显著性水平(通常为 0.05),则可以拒绝零假设(即系数等于零),认为该系数在统计上是显著的。在 statsmodels 的回归结果中,pvalue 属性包含了每个回归系数的 p 值。

从模型I可以看到,ACT变量的pvalue为0.01<0.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显著影响大学成绩GPA,params为0.45代表每1单位的高考成绩的增加将导致大学GPA平均增加0.45个单位。

(2)分类分析:

以ST——我国股市特有的一项保护投资者利于的决策为例,当上市公司因财务状况不佳导致投资者难以预测其前景时,交易所会标记该公司股票为ST,并进行一系列限制措施。我们想研究被ST的公司其背后的因素,并尝试通过利用公司的财务指标提前预测某上市公司在未来是否会被ST,是一件很有意义的举措。而在这项任务中,因变量就是公司是否会被ST。该例中自变量是一些财务指标,如ARA、ASSET等。

在这段代码中,使用了 statsmodels 库来进行逻辑回归分析。

  1. st_logit=sm.formula.logit('ST~ARA+ASSET+ATO+ROA+GROWTH+LEV+SHARE',data=ST).fit()

    • sm.formula.logit 是 statsmodels 中用于逻辑回归(Logistic Regression)的函数。
    • 'ST~ARA+ASSET+ATO+ROA+GROWTH+LEV+SHARE' 是一个公式字符串,表示模型的因变量是 ST,自变量包括 ARAASSETATOROAGROWTHLEV 和 SHARE。这里的 ~ 符号用于分隔因变量和自变量。
    • data=ST 指定了数据源,其中 ST 应该是一个包含相应变量的 pandas DataFrame。
    • .fit() 是用来拟合模型的方法。
  2. print(st_logit.summary())

st_logit.summary() 会输出一个格式化的表格,包含以下主要部分:

  • Model: 模型类型和方法,这里显示为逻辑回归。
  • Method: 拟合方法,通常是最大似然估计。
  • No. Observations: 观测值的数量。
  • Df Residuals: 自由度,即观测值减去参数数量。
  • Df Model: 模型自由度,即自变量的数量。
  • Pseudo R-squared: 伪 R-squared 值,用于衡量模型的拟合优度,但由于逻辑回归不是线性回归,这个值与普通线性回归的 R-squared 值有所不同。
  • Log-Likelihood: 对数似然函数值,用于衡量模型的拟合优度。
  • AIC: 赤池信息准则,用于模型选择,越小表示模型越好。
  • BIC: 贝叶斯信息准则,也是用于模型选择,考虑了模型复杂度的影响。
  • LL-Null: 零模型的对数似然函数值,即只包含截距项的模型。
  • LLR p-value: 似然比检验的 p 值,用于检验模型中所有自变量的联合显著性。

参数估计值

  • coef: 参数的估计值,即每个自变量的系数。
  • std err: 标准误差,表示参数估计值的不确定性。
  • z: z 值,标准误差的系数估计值的比值。
  • P>|z|: p 值,用于检验每个参数是否显著不同于零。

置信区间

  • [0.025] 和 [0.975] 分别表示 95% 置信区间的下限和上限。
# 加载基础包
import pandas as pd
import numpy as np
import statsmodels.api as sm
from scipy import stats

# 读取数据
ST=pd.read_csv('./data/ST.csv')
ST.head()

st_logit=sm.formula.logit('ST~ARA+ASSET+ATO+ROA+GROWTH+LEV+SHARE',data=ST).fit()
print(st_logit.summary())

样例输出结果分析:

 其中pvalue为第四列标记为P(>| z |)的数值。与回归分析类似,一般认为,p值小于0.05的自变量是显著的,也就是说,有统计学证据表明该变量会影响因变量为1的概率。

注:

在逻辑回归中,"P>|z|" 表示的是每个回归系数的 p 值,它用于检验该系数是否在统计上显著地不同于零。具体来说,p 值是观察到的统计量(在逻辑回归中通常是 Wald 统计量)或更极端情况下的概率。以下是如何解释 "P>|z|" 的结果:

  1. p 值的定义

    • p 值是当零假设(即系数等于零)为真时,观察到的统计量或更极端情况出现的概率。在逻辑回归中,零假设是某个自变量对因变量没有影响。
  2. 显著性水平

    • 通常,显著性水平(α)被设定为 0.05。这意味着如果 p 值小于 0.05,我们有足够的证据拒绝零假设,认为该系数在统计上显著不同于零。
  3. 解释 "P>|z|"

    • 如果 "P>|z|" 的值小于显著性水平(例如 0.05),则可以认为该系数显著,即该自变量对因变量有显著影响。
    • 如果 "P>|z|" 的值大于显著性水平(例如 0.05),则没有足够的证据拒绝零假设,认为该系数不显著,即该自变量对因变量没有显著影响。
  4. z 值

    • z 值是系数估计值与标准误差的比值。它是一个标准化的统计量,用于计算 p 值。z 值越大,p 值越小,系数显著的可能性越大。
  5. 逻辑回归中的 Wald 检验

    • 在逻辑回归中,Wald 检验是一种常用的方法来检验单个系数的显著性。它基于系数估计值的标准误差和系数的 z 值来计算 p 值。
  6. 模型整体的显著性

    • 除了检验单个系数的显著性外,还可以通过似然比检验(Likelihood Ratio Test, LRT)来检验模型整体的显著性。这通常通过比较模型的对数似然函数值(LL)和零模型的对数似然函数值(LL-Null)来实现。
3.2 假设检验

本质:根据样本信息与已知信息,对一个描述总体性质的命题进行“是或否”的检验与回答,假设检验验证的不是样本本身的性质,而是样本所在总体的性质

假设检验大体上可分为两种:参数假设检验与非参数假设检验。若假设是关于总体的一个参数或者参数的集合,则该假设检验就是参数假设检验,刚刚的例子的假设是关于总体均值的,均值是参数,因此这是一个参数假设检验;若假设不能用一个参数集合来表示,则该假设检验就是非参数假设检验,一个典型就是正态性检验。

3.2.1 正态性检验

参数检验比非参数检验更灵敏,因此一旦数据是正态分布的,我们应该使用参数检验,此时对数据进行正态性检验就非常有必要了。

(1)判断数据的正态性的三种方法

①可视化判断-正态分布概率图

# 生成1000个服从正态分布的数据
data_norm = stats.norm.rvs(loc=10, scale=10, size=1000) # rvs(loc,scale,size):生成服从指定分布的随机数,loc:期望;scale:标准差;size:数据个数
# 生成1000个服从卡方分布的数据
data_chi=stats.chi2.rvs(2,3,size=1000)

# 画出两个概率图
fig=plt.figure(figsize=(12,6))
ax1=fig.add_subplot(1,2,1)
plot1=stats.probplot(data_norm,plot=ax1) # 正态数据
ax2=fig.add_subplot(1,2,2)
plot2=stats.probplot(data_chi,plot=ax2) # 卡方分布数据

其中,

  • stats.norm.rvs 是 scipy.stats 库中用于生成正态分布随机数的函数
  • loc=10:正态分布的均值(期望值)为 10。
  • scale=10:正态分布的标准差为 10。
  • size=1000:生成 1000 个随机数。
  • stats.chi2.rvs 是 scipy.stats 库中用于生成卡方分布随机数的函数。
  • 2:卡方分布的自由度为 2。
  • 3:卡方分布的缩放因子,通常用于生成非标准卡方分布。
  • size=1000:生成 1000 个随机数。
  • fig=plt.figure(figsize=(12,6)):创建一个图形窗口,大小为 12x6 英寸。
  • ax1=fig.add_subplot(1,2,1):在图形窗口中添加一个子图,布局为 1 行 2 列,这是第一个子图。
  • plot1=stats.probplot(data_norm,plot=ax1):在第一个子图中绘制正态分布数据的概率图。
  • ax2=fig.add_subplot(1,2,2):在图形窗口中添加第二个子图。
  • plot2=stats.probplot(data_chi,plot=ax2):在第二个子图中绘制卡方分布数据的概率图。

样例效果分析:

  • 正态分布、卡方分布的概率图

    • 概率图显示了数据的分位数与正态分布/卡方分布的理论分位数之间的关系。
    • 如果数据确实服从正态分布/卡方分布,那么图中的点应该大致沿着一条直线排列。
  • 从图中可以看到,正态性数据在正态分布概率图中十分接近直线y=x,而卡方分布数据则几乎完全偏离了该直线。可见,概率图可视化确实可以起到初步判断数据正态性的作用。实际应用中,由于数据的复杂性,仅使用一种方法判断正态性有可能产生一定的误差,因此我们通常同时使用多种方法进行判断。如果不同方法得出的结论不同,此时就需要仔细观察数据的特征,寻找结果不一致的原因。如:若Shapiro-Wilk test显著(非正态),D'Agostino's K-squared test不显著(正态),则有可能是因为样本量较大,或者样本中存在重复值现象,如果事实确实如此,那么我们就应该采纳D'Agostino's K-squared test的结论而非Shapiro-Wilk test的结论。(在这里,我们假设数据服从正态分布,p<0.05则拒绝假设,换个说法:p<0.05数据不服从正态分布,如果p>=0.05说明没有证据说明数据不服从正态分布)
  • 注:
  • Shapiro-Wilk检验:适用于小样本数据,当怀疑数据不服从正态分布时,可以使用此检验。
  • D'Agostino's K-squared检验:适用于大样本数据,当需要同时考虑数据的偏度和峰度时,可以使用此检验。

②Shapiro-Wilk检验

  1. 定义:Shapiro-Wilk检验是一种统计检验,用于检验一组数据是否服从正态分布。该检验是由 Samuel S. Shapiro 和 Martin Wilk 在1965年提出的。

  2. 原理:检验的基本思想是将数据的顺序统计量与正态分布的顺序统计量进行比较。如果数据来自正态分布,那么这些顺序统计量应该与正态分布的理论值接近。

  3. 计算:检验统计量 W 是基于数据的均值和标准差计算的,具体计算公式较为复杂,通常通过统计软件实现。

  4. p 值:检验结果会给出 p 值。如果 p 值小于显著性水平(通常为 0.05),则拒绝正态分布的零假设,认为数据不服从正态分布。

  5. 优点:Shapiro-Wilk检验对于小样本(n<50)非常有效,并且对正态分布的偏离非常敏感。

  6. 缺点:随着样本量的增加,检验的功效会降低。

③D'Agostino's K-squared检验

  1. 定义:D'Agostino's K-squared检验是由 Ralph B. D'Agostino 在1971年提出的,用于检验数据是否服从正态分布。

  2. 原理:检验基于数据的偏度(Skewness)和峰度(Kurtosis)。偏度描述了数据分布的对称性,峰度描述了数据分布的尖峭程度。注:

    偏度(Skewness)

    峰度(Kurtosis)
    定义描述数据分布不对称性的统计量。它衡量数据分布相对于正态分布的偏斜程度。描述数据分布峰态的统计量,衡量数据分布的尖峭程度或尾部的厚度。
    计算n 是样本大小,xi​ 是每个观测值,xˉ 是样本均值,s 是样本标准差。n 是样本大小,xi​ 是每个观测值,xˉ 是样本均值,s 是样本标准差
    解释

    正偏度:如果偏度大于 0,数据分布的尾部向右延伸,称为正偏(或右偏);

    负偏度:如果偏度小于 0,数据分布的尾部向左延伸,称为负偏(或左偏);

    无偏度:如果偏度接近 0,数据分布接近对称

    高峰度:如果峰度大于 3(或大于正态分布的峰度),数据分布比正态分布更尖峭,尾部更薄,称为高峰度(或尖峰);

    低峰度:如果峰度小于 3,数据分布比正态分布更平坦,尾部更厚,称为低峰度(或扁峰);

    正态峰度:如果峰度接近 3,数据分布的峰态接近正态分布。

    意义帮助识别数据的对称性问题,如收入分布的正偏帮助识别数据的极端值,如金融市场的极端波动
  3. 优点:D'Agostino's K-squared检验对大样本(n>50)更为有效,并且可以提供关于数据分布形状的额外信息。

  4. 缺点:对于小样本,检验的功效不如 Shapiro-Wilk检验。

3.2.2 单组样本均值假定的检验
①t检验——用于判断是否拒绝零假设
  1. 定义:t检验是一种用于比较两组数据均值是否存在显著差异的统计方法,通常用于样本量较小且数据近似正态分布的情况。

  2. 类型:单样本t检验用于比较单个样本的均值与已知的总体均值是否存在显著差异;独立样本t检验(双样本t检验)用于比较两组独立样本的均值是否存在显著差异;配对样本t检验用于比较同一组样本在不同条件下的均值是否存在显著差异。

  3. 原理:t检验基于t分布,通过计算t统计量来评估样本均值与假设值之间的差异是否显著。t统计量的计算公式为:其中,xˉ 是样本均值,μ0​ 是假设的总体均值,s 是样本标准差,n 是样本大小。

  4. p值:检验结果会给出p值,用于判断是否拒绝零假设(即样本均值与假设值无显著差异)。如果p值小于显著性水平(通常为0.05),则拒绝零假设。

  5. 应用:适用于样本量较小且数据近似正态分布的情况,是社会科学、生物科学等领域中常用的统计方法。

代码示例:

pVal=pd.Series(dtype='float64')
  • 使用pandas.Series初始化一个空的序列,用于存储p值。
_, pVal['t-test'] = stats.ttest_1samp(data, checkvalue,alternative=alternative)
  • 使用scipy.stats.ttest_1samp函数进行单样本t检验。该函数比较样本均值与checkvalue
  • 第一个返回值(通常是t统计量)被丢弃(用下划线_表示)。
  • 第二个返回值是p值,存储在pVal['t-test']中。

 ②Wilcoxon符号秩检验——用于判断是否拒绝零假设

  1. 定义:Wilcoxon符号秩检验是一种非参数检验方法,用于比较两个独立样本的中位数是否存在显著差异,或比较单个样本的中位数与已知中位数是否存在显著差异。

  2. 原理:该检验基于秩次(即数据在排序后的位置),而不是数据的具体数值。它适用于不满足正态分布假设的数据,或样本量较小的情况。

  3. 计算:首先计算每个观测值与零假设值(或另一组样本的中位数)的差值,然后对差值的绝对值进行排序,并计算正负差值的秩和。

  4. p值:检验结果同样会给出p值,用于判断是否拒绝零假设(即两组样本的中位数无显著差异,或样本中位数与假设中位数无显著差异)。如果p值小于显著性水平,则拒绝零假设。

  5. 应用:适用于数据不满足正态分布假设的情况,或样本量较小且数据类型为序数数据的情况。常用于医学研究、心理学研究等领域。

_, pVal['wilcoxon'] = stats.wilcoxon(data-checkvalue,alternative=alternative)
  • 使用scipy.stats.wilcoxon函数进行Wilcoxon符号秩检验。该函数比较样本数据与checkvalue的差值。
  • 同样,第一个返回值(通常是Wilcoxon秩和)被丢弃。
  • 第二个返回值是p值,存储在pVal['wilcoxon']中。
3.2.3 两组样本的均值相等性检验

前要工作:在进行两组之间的均值比较之前,需要进行重要的判断:这两组样本之间是否独立。这个问题的答案将决定着使用何种类型的检验。事实上,两个样本是否在抽样意义上独立,是没有固定答案的,因为在很多情况下,我们既不能保证两个样本间完全独立,也很难判断出两者是否存在相关性。我们只能在抽样的时候使用更加科学的抽样方法,尽可能地避免样本间的相互影响。若两个样本的总体都服从正态分布,那么我们可以使用双样本t检验。如果不服从正态分布,则可以使用Mann-whitney u秩和检验,Mann-whitney u秩和检验是一种非参数检验。

双样本t检验(独立样本t检验)

Mann-whitney u秩和检验
定义一种用于比较两组独立样本均值是否存在显著差异的统计方法一种非参数检验方法,用于比较两组独立样本的分布是否存在显著差异
适用条件
  • 两组样本是独立的。
  • 样本数据近似正态分布,或者样本量足够大(通常认为大于30)时,根据中心极限定理,样本均值的分布会接近正态分布。
  • 两组样本是独立的。
  • 不需要满足正态分布假设,适用于数据分布未知或非正态分布的情况
计算其中,xˉ1​ 和 xˉ2​ 是两组样本的均值,μ1​ 和 μ2​ 是两组样本的总体均值(通常假设相等),s1​ 和 s2​ 是两组样本的标准差,n1​ 和 n2​ 是两组样本的大小。
  • 将两组数据合并并进行排序,赋予秩次(相同的值赋予平均秩次)。
  • 计算每组数据的秩和(即每组数据中所有观测值的秩次之和)。、
  • 注:
  • 秩次是一个数据在一组数据中的相对位置或顺序。例如,在一组按大小排序的数据中,最大的数据的秩次是1,次大的秩次是2,依此类推。
  • 如果数据中有相同的值,则这些值会共享一个秩次,通常取它们秩次的平均值。例如,如果有3个相同的值,它们的秩次可以是 (3 + 4 + 5) / 3 = 4。
  • 秩和是一组数据中所有数据点秩次的总和。对于一组数据,首先计算每个数据点的秩次,然后将这些秩次相加得到秩和。
  • n1​ 和 n2​ 是两组样本的大小,R是较小样本的秩和。
  • 根据U值和样本大小,查找Mann-Whitney U检验表或使用统计软件计算p值。
p值
  • 根据t统计量和自由度(n1​+n2​−2),查找t分布表或使用统计软件计算p值。
  • 如果p值小于显著性水平(通常为0.05),则拒绝零假设,认为两组样本均值存在显著差异。
  • 根据秩和和样本大小,查找Mann-Whitney U检验表或使用统计软件计算p值。
  • 如果p值小于显著性水平,则拒绝零假设,认为两组样本的分布存在显著差异。
3.2.4 方差分析-多组样本间的均值相等性检验

方差分析通过分解数据中的总变异(Total Variation)为由不同因素引起的变异和随机误差两部分,来检验不同组别间均值的差异。总变异可以表示为:

总变异=组间变异+组内变异+误差总变异=组间变异+组内变异+误差

  1. 组间变异(Between Groups Variation):不同组别间的差异引起的变异。
  2. 组内变异(Within Groups Variation):同一组内个体间的差异引起的变异。
  3. 误差:未被解释的随机变异。
方差分析的类型
  1. 单因素方差分析(One-Way ANOVA):只涉及一个因素(自变量),该因素有多个水平(组别)。用于检验不同组别间均值是否存在显著差异。

  2. 双因素方差分析(Two-Way ANOVA):涉及两个因素(自变量),每个因素有多个水平。用于检验两个因素及其交互作用对结果变量的影响。

  3. 多因素方差分析(Factorial ANOVA):涉及两个以上的自变量,每个自变量有多个水平。可以检验多个因素及其交互作用对结果变量的影响。

  4. 重复测量方差分析(Repeated Measures ANOVA):用于同一样本在不同时间点或不同条件下的测量数据。可以检验时间效应、条件效应及其交互作用。

方差分析的步骤

  1. 建立假设:零假设(H0):所有组别的均值相等;备择假设(H1):至少有两个组别的均值不相等。

  2. 计算组间和组内平方和:组间平方和(SSB):衡量组间变异的大小;组内平方和(SSW):衡量组内变异的大小。

  3. 计算均方误差:组间均方误差(MSB):MSB=SSB/dfB;组内均方误差(MSW):MSW=SSW/dfW,其中,dfB和dfW分别是组间和组内的自由度。

  4. 计算F统计量: 𝐹=MSB/MSW F统计量用于衡量组间变异与组内变异的比率。

  5. 确定显著性水平并计算p值

    • 根据F统计量和相应的自由度,查找F分布表或使用统计软件计算p值。
    • 如果p值小于显著性水平(通常为0.05),则拒绝零假设,认为至少有两个组别的均值存在显著差异。
  6. 多重比较如果方差分析显示显著差异,可以进一步进行多重比较,如Tukey、Bonferroni或Scheffé方法,来确定哪些组别之间的差异显著。

        尽管在大样本下,非正态性数据的方差分析也是稳健的,但是在小样本下,对非正态性数据做方差分析还是可能存在误差。此时,我们可以使用kruskalwallis检验。

Part 4 随机过程与随机模拟

随机过程是概率论在时间或空间维度上的扩展,它提供了一种研究随时间或空间变化的随机现象的框架。通过随机过程,可以更深入地理解复杂系统中的随机性和不确定性。

4.1 随机过程的定义
  1. 定义:随机过程是定义在概率空间上的随机变量序列,通常表示为 𝑋(𝑡),其中 𝑡可以是时间或空间的参数。每个 𝑡对应一个随机变量 𝑋(𝑡),整个集合 {𝑋(𝑡),𝑡∈𝑇}构成了一个随机过程,其中 𝑇是参数集,可以是实数集或其他集合。

  2. 特点随机性:每个随机变量 𝑋(𝑡) 都具有随机性,即在相同的 𝑡下,𝑋(𝑡)可以取不同的值;时间或空间依赖性:随机过程的值不仅依赖于 𝑡,还可能依赖于过程的历史,即 𝑋(𝑡) 的值可能与 𝑋(𝑠)(𝑠<𝑡)的值有关。

  3. 分类离散时间随机过程:参数 𝑡 取离散值,如整数集合;连续时间随机过程:参数 𝑡取连续值,如实数集合。

  4. 重要性质

  • 概率分布:研究随机过程在不同时间点的分布,如 𝑋(𝑡)的概率密度函数或分布函数。
  • 联合分布:研究随机过程在多个时间点的联合分布,如 (𝑋(𝑡1),𝑋(𝑡2),…,𝑋(𝑡𝑛))的联合分布。
  • 期望和方差:研究随机过程的期望值和方差,以及它们随时间的变化。
4.2 概率论与随机过程的区别
随机过程概率论
研究对象研究随时间或空间变化的随机现象,关注的是一系列随机变量的集合及其统计规律。主要研究单个随机事件的概率和随机变量的概率分布。它关注的是单个随机变量的性质,如期望、方差等。
时空因素考虑时间或空间的变化,研究随机变量随时间或空间的演变。不考虑时间或空间的变化,主要关注静态的随机变量。
数学工具除了使用概率论的工具外,还引入了如协方差、自相关函数、谱密度等概念来描述随机过程的性质。使用概率、分布函数、期望、方差等基本概念和工具。
应用场景应用于随时间或空间变化的随机现象的分析,如股票价格的波动、信号的传输等。应用于单个随机事件的分析,如掷骰子、抽卡片等。

        随机过程是概率论在时间或空间维度上的扩展,它提供了一种研究随时间或空间变化的随机现象的框架。通过随机过程,可以更深入地理解复杂系统中的随机性和不确定性。

Part 5 数据可视化

5.1 什么是一个好的数据可视化图表?
  • 图表展示的信息全面且无歧义
  • 图表表达的信息越多、越全面越好
  • 通俗易懂,不能太专业
5.2 python数据可视化工具简介
5.2.1 Matplotlib

是一个面向对象的绘图工具库,pyplot是Matplotlib最常用的一个绘图接口,提供了一整套和Matlab相似的命令API,适合交互式制图,还可以将它作为绘图控件,嵌入其它应用程序中。

import matplotlib.pyplot as plt
  • 创建一个图形对象,并设置图形对象的大小(可以想象成在白纸中添加一个图,并设置图的大小):plt.figure(figsize=(6,4))
  • 在纸上的坐标系中绘制散点:plt.scatter(x=x, y=y)
  • 设置x轴的标签label:plt.xlabel('x')
  • 设置y轴标签的label:plt.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(x=x, y=y)
plt.xlabel('x')
plt.ylabel('y')
plt.title('y = sin(x)')
plt.show()

 使用一个例子体会Matplotlib的面向对象绘图的特性:

【例子】绘制y = sin(x) 和 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), num='first')
fig1.suptitle('y = sin(x)')
plt.scatter(x=x, y=y1)
plt.xlabel('x')
plt.ylabel('y')

# 绘制第二个图:
fig2 = plt.figure(figsize=(6,4), num='second')
fig2.suptitle('y = cos(x)')
plt.scatter(x=x, y=y2)
plt.xlabel('x')
plt.ylabel('y')

plt.show()

5.2.2 Seaborn

Seaborn主要用于统计分析绘图的,它是基于Matplotlib进行了更高级的API封装。Seaborn比matplotlib更加易用,尤其在统计图表的绘制上,因为它避免了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",data=df)

Seaborn把数据拟合等统计属性高度集成在绘图函数中,绘图的功能还是构筑在Matplotlib之上。因此,Seaborn相当于是完善了统计图表的Matplotlib工具库,二者应该是相辅相成的。因此,在实际的可视化中,往往一起使用Matplotlib和Seaborn。

5.2.3 Plotnine

Plotnine 是一个基于 Python 的图形库,它提供了一种类似于 R 语言中 ggplot2 包的图形语法,旨在为数据科学家和分析师提供一个简单、灵活且功能强大的工具来创建高质量的图形和图表。

5.3 基本图标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)'], width=0.6, align='center', orientation='vertical', 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",水平条:"horizontal"
"""
plt.title("不同城市的房价对比图")   # 在axes1设置标题
plt.xlabel("城市")    # 在axes1中设置x标签
plt.ylabel("房价/w")    # 在axes1中设置y标签
plt.grid(b=True, which='both')  # 在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)'],fontsize=13)   # 添加数据注释
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(x='city', y='house_price(w)', fill='city', label='house_price(w)'))+
    geom_bar(stat='identity')+
    labs(x="城市", y="房价(w)", title="不同城市的房价对比图")+
    geom_text(nudge_y=0.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)'], 
    width=0.3, 
    align='center', 
    orientation='vertical', 
    label='年份:2021'
    )
plt.bar(
    np.arange(len(np.unique(data['城市'])))+0.15, 
    data.loc[data['年份']==2022,'房价(w)'], 
    width=0.3, 
    align='center', 
    orientation='vertical', 
    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(b=True, which='both')  # 在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],fontsize=13)   # 添加数据注释

data_2022 = data.loc[data['年份']==2022,:]
for i in range(len(data_2022)):
    plt.text(i+0.15-0.05, data_2022.iloc[i,2]+0.05, data_2022.iloc[i,2],fontsize=13)   # 添加数据注释
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['年份'], ordered=True, categories=data['年份'].unique())
p_mult_bar = (
    ggplot(data, aes(x='城市', y='房价(w)', fill='年份'))+
    geom_bar(stat='identity',width=0.6, position='dodge')+
    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]
})
tmp=data.set_index(['城市','年份'])['房价(w)'].unstack()
data=tmp.rename_axis(columns=None).reset_index()
data.columns = ['城市','2021房价','2022房价']
print(data)

plt.figure(figsize=(10,6))
plt.bar(
    data['城市'], 
    data['2021房价'], 
    width=0.6, 
    align='center', 
    orientation='vertical', 
    label='年份:2021'
    )
plt.bar(
    data['城市'], 
    data['2022房价'], 
    width=0.6, 
    align='center', 
    orientation='vertical', 
    bottom=data['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]
})
tmp=data.set_index(['城市','年份'])['房价(w)'].unstack()
data=tmp.rename_axis(columns=None).reset_index()
data.columns = ['城市','2021房价','2022房价']
print(data)

plt.figure(figsize=(10,6))
plt.bar(
    data['城市'], 
    data['2021房价']/(data['2021房价']+data['2022房价']), 
    width=0.4, 
    align='center', 
    orientation='vertical', 
    label='年份:2021'
    )
plt.bar(
    data['城市'], 
    data['2022房价']/(data['2021房价']+data['2022房价']), 
    width=0.4, 
    align='center', 
    orientation='vertical', 
    bottom=data['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), endpoint=False)   # 每个坐标点的位置
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-', lw=2, label='艾希') # 画出雷达图的点和线
plt.fill(theta, aixi, facecolor='r', alpha=0.5) # 填充
plt.polar(theta, nuoshou, 'bo-', lw=2, label='诺手')  # 画出雷达图的点和线
plt.fill(theta, nuoshou, facecolor='b', alpha=0.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, c='dodgerblue', marker=".", s=50)
plt.xlabel("x")
plt.ylabel("y1")
plt.title("y1与x不存在关联关系")

plt.subplot(2,2,2)  #创建两行两列的子图,并绘制第二个子图
plt.scatter(x, y2, c='tomato', marker="o", s=10)
plt.xlabel("x")
plt.ylabel("y2")
plt.title("y2与x存在关联关系")

plt.subplot(2,2,3)  #创建两行两列的子图,并绘制第三个子图
plt.scatter(x, y3, c='magenta', marker="o", s=10)
plt.xlabel("x")
plt.ylabel("y3")
plt.title("y3与x存在关联关系")

plt.subplot(2,2,4)  #创建两行两列的子图,并绘制第四个子图
plt.scatter(x, y4, c='deeppink', marker="s", s=10)
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(x='x', y='y', fill='class', shape='class'), color='black', size=2)+
    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 sm

from sklearn.preprocessing import PolynomialFeatures  # 构造多项式
x = np.linspace(-10, 10, 100)
y = np.square(x) + np.random.randn(100)*100
x_poly2 = PolynomialFeatures(degree=2).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, frac=2/3)[:, 1]

plt.figure(figsize=(8, 8))
plt.subplot(2,2,1)
plt.scatter(x, y, c='tomato', marker="o", s=10)
plt.plot(x, y_linear_pred, c='dodgerblue')
plt.xlabel("x")
plt.ylabel("y")
plt.title("带线性趋势线的散点图")

plt.subplot(2,2,2)
plt.scatter(x, y, c='tomato', marker="o", s=10)
plt.plot(x, y_poly_pred, c='dodgerblue')
plt.xlabel("x")
plt.ylabel("y")
plt.title("带二次趋势线的散点图")

plt.subplot(2,2,3)
plt.scatter(x, y, c='tomato', marker="o", s=10)
plt.plot(x, y_exp_pred, c='dodgerblue')
plt.xlabel("x")
plt.ylabel("y")
plt.title("带指数趋势线的散点图")

plt.subplot(2,2,4)
plt.scatter(x, y, c='tomato', marker="o", s=10)
plt.plot(x, y_loess_pred, c='dodgerblue')
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, columns=feature)
df['label'] = label

label_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(x=df_label['sepal length (cm)'], y=df_label['sepal width (cm)'], s=20, label=i)
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_vars='index', var_name='variable', value_name='corr_xy') #将矩阵宽数据变成长数据
df['abs_corr'] = np.abs(df['corr_xy'])
p1 = (
    ggplot(df, aes(x='index', y='variable', fill='corr_xy', size='abs_corr'))+
    geom_point(shape='o', color='black')+
    scale_size_area(max_size=11, guide=False)+
    scale_fill_cmap(name='RdYIBu_r')+
    coord_equal()+
    labs(x="Variable", y="Variable")+
    theme(dpi=100, figure_size=(4.5,4.55))
)
p2 = (
    ggplot(df, aes(x='index', y='variable', fill='corr_xy', size='abs_corr'))+
    geom_point(shape='s', color='black')+
    scale_size_area(max_size=10, guide=False)+
    scale_fill_cmap(name='RdYIBu_r')+
    coord_equal()+
    labs(x="Variable", y="Variable")+
    theme(dpi=100, figure_size=(4.5,4.55))
)
p3 = (
    ggplot(df, aes(x='index', y='variable', fill='corr_xy', label='corr_xy'))+
    geom_tile(color='black')+
    geom_text(size=8, color='white')+
    scale_fill_cmap(name='RdYIBu_r')+
    coord_equal()+
    labs(x="Variable", y="Variable")+
    theme(dpi=100, 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'], bins=20, alpha=0.85)
plt.xlabel("mpg")
plt.ylabel("count")
plt.show()

# 使用matplotlib绘制箱线图
import seaborn as sns 
from plotnine.data import mtcars

data = mtcars
data['carb'] = data['carb'].astype('category')
plt.figure(figsize=(8, 6))
sns.boxenplot(x='carb', y='mpg', data=mtcars, linewidth=0.2, palette=sns.husl_palette(6, s=0.9, l=0.65, h=0.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='己方输出', ascending=False).reset_index()
df_other = df[['敌方', '敌方输出']].sort_values(by='敌方输出', ascending=False).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, startangle=90, shadow=True, colors=color_list, labels=df_our['己方'].tolist(), explode=(0,0,0,0,0.3), autopct='%.2f%%')


plt.subplot(1,2,2)
plt.pie(df_other['敌方输出'].values, startangle=90, shadow=True, colors=color_list, labels=df_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='己方输出', ascending=False).reset_index()
df_other = df[['敌方', '敌方输出']].sort_values(by='敌方输出', ascending=False).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, startangle=90, shadow=True, colors=color_list, wedgeprops=wedgeprops, labels=df_our['己方'].tolist(), explode=(0,0,0,0,0.3), autopct='%.2f%%')
plt.text(0, 0, '己方' , ha='center', va='center', fontsize=30)

plt.subplot(1,2,2)
plt.pie(df_other['敌方输出'].values, startangle=90, shadow=True, colors=color_list, wedgeprops=wedgeprops, labels=df_other['敌方'].tolist(), explode=(0,0,0,0,0.3), autopct='%.2f%%')
plt.text(0, 0, '敌方' , ha='center', va='center', fontsize=30)
plt.show()

 (4)时间序列型图表
# 使用Plotnine绘制时间序列线图:
df = pd.read_csv(
    './data/AirPassengers.csv'
    )  # 航空数据1949-1960
df['date'] = pd.to_datetime(df['date'])
p1 = (
    ggplot(df, aes(x='date', y='value'))+
    geom_line(size=1, color='red')+
    scale_x_date(date_labels="%Y", date_breaks="1 year")+
    xlab('date')+
    ylab('value')
)
print(p1)

# 使用Matplotlib绘制时间序列折线图
plt.figure(figsize=(8,6))
plt.plot(df['date'], df['value'], color='red')
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(x='date_list', y='value_list', group='factor(Class)', color='factor(Class)'))+
    geom_line(size=1)+
    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'], color='red', alpha=0.86, label='value1')
plt.plot(data['date_list'], data['value_list2'], color='blue', alpha=0.86, label='value2')
plt.legend()
plt.xlabel('date')
plt.ylabel('value')
plt.show()

Part 6 插值模型

插值方法实际上除了是用于处理缺失值一个很好的方法,也是做数据填充的很好的方法。

  1. 语法一致性:Plotnine 采用了类似于 ggplot2 的图形语法,使得从 R 语言迁移到 Python 的用户可以轻松上手。通过使用 + 操作符,可以逐步构建图形,添加图层、坐标轴、标签等。

  2. 数据驱动:Plotnine 强调数据驱动的图形设计,允许用户通过数据框(如 pandas DataFrame)来定义图形的各个方面。

  3. 灵活性:支持多种图形类型,包括散点图、线图、条形图、箱线图、直方图等;允许用户自定义图形的各个方面,如颜色、形状、大小等。

  4. 扩展性:可以通过安装额外的包来扩展 Plotnine 的功能,例如使用 plotnine.facet 进行小面分割,或使用 plotnine.scales 进行更复杂的数据缩放。

  5. 交互性:虽然 Plotnine 本身不支持交互式图形,但可以通过与 Jupyter Notebook 结合使用或其他工具来实现交互式图形的展示。

import pandas as pd
from plotnine import ggplot, aes, geom_line, geom_point
from plotnine import *     # 讲Plotnine所有模块引入
from plotnine.data import *   # 引入PLotnine自带数据集
mpg.head()


# 创建数据框
data = pd.DataFrame({
    'x': [1, 2, 3, 4],
    'y': [1, 4, 9, 16]
})

# 创建图形
p = (ggplot(data, aes(x='x', y='y'))
     + geom_line() 
     + geom_point()
     + ggtitle("Simple Plotnine Plot"))

# 显示图形
print(p)

# 绘制汽车在不同驱动系统下,发动机排量与耗油量的关系
p1 = (
    ggplot(mpg, aes(x='displ', y='hwy', color='drv'))     # 设置数据映射图层,数据集使用mpg,x数据使用mpg['displ'],y数据使用mpg['hwy'],颜色映射使用mog['drv']
    + geom_point()       # 绘制散点图图层
    + geom_smooth(method='lm')        # 绘制平滑线图层
    + labs(x='displacement', y='horsepower')     # 绘制x、y标签图层
)
print(p1)   # 展示p1图像
6.1 线性插值法

线性插值对两个点中的解析式是按照线性方程来建模。比如我们得到的原始数据列{y}和数据的下标{x},这里数据下标x可能并不是固定频率的连续取值而是和y一样存在缺失的。给定了数据点(xk,yk)和(xk+1,yk+1),需要对两个点之间构造直线进行填充。

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, k=1)  # 样本点导入,生成参数
iy1 = spi.splev(new_x, ipo1)  # 根据观测点和样条参数,生成插值
  • 这段代码的主要目的是:

  • 生成一组正弦函数的样本点。
  • 定义一组更密集的插值点。
  • 使用一阶样条插值方法,在插值点上计算插值结果。
  • 其中:
  • numpy库用于数值计算,提供了大量的数学函数和矩阵操作。
  • scipy.interpolate库用于插值,提供了各种插值方法。
  • np.arange函数生成一个从 -2 * np.pi 到 2 * np.pi 的等差数列,步长为1。这里生成了从 -2π 到  的整数倍的π值。
  • np.sin函数计算每个样本点X的正弦值,生成样本点Y。
  • 定义了一个新的数组 new_x,用于插值的点。与样本点X不同,这里的步长为0.1,即在 -2π 到  之间每隔0.1生成一个点。
  • spi.splrep函数用于生成样条插值的参数。参数 k=1 表示使用一阶样条(线性样条)。该函数接受三个参数:
    • X:样本点的横坐标。
    • Y:样本点的纵坐标。
    • k:样条的阶数,这里使用一阶样条。
  • spi.splev函数用于根据生成的样条参数和新的插值点 new_x 进行插值。该函数返回插值后的纵坐标值。
 6.2 三次样条插值

三次样条插值是一种非常自然的插值方法。它将两个数据点之间的填充模式设置为三次多项式。

它假设在数据点(xk,yk)和(xk+1,yk+1)之间的三次式叫做Ik,通过解方程的形式可以解出每一只三次式。简而言之,某个数据点前后两条三次函数不仅在当前的函数值相等,一次导数值和二次导数值也要保持相等。

# 进行三次样条拟合
ipo3 = spi.splrep(X, Y, k=3)  # 样本点导入,生成参数
iy3 = spi.splev(new_x, ipo3)  # 根据观测点和样条参数,生成插值
6.3 拉格朗日插值

def lagrange(x0,y0,x):
    y=[]
    for k in range(len(x)):
        s=0
        for i in range(len(y0)):
           t=y0[i]
           for j in range(len(y0)):
              if i!=j:
                t*=(x[k]-x0[j])/(x0[i]-x0[j])
           s+=t
        y.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]
x=np.arange(1,5,1/30)
f1=interp1d(x0,y0,'linear')
y1=f1(x)
f2=interp1d(x0,y0,'cubic')
y2=f2(x)
y3=lagrange(x0,y0,x)
plt.plot(x0,y0,'r*')
plt.plot(x,y1,'b-',x,y2,'y-',x,y3,'r-')
plt.show()

 6.4 三种插值方法的详细介绍和比较

①线性插值

  • 简单:计算简单,易于实现。
  • 局限性:只能保证在已知点处的插值准确,中间点的插值效果可能不理想,尤其是在数据变化剧烈的情况下。

②三次样条插值

  • 平滑性:插值曲线在数据点处不仅值连续,一阶导数也连续,因此曲线平滑。
  • 局部性:每个多项式仅依赖于附近的数据点,对远离的数据点影响较小。
  • 计算复杂性:相对于线性插值,计算更为复杂,需要解决一个线性方程组。

③拉格朗日插值

  • 精确性:对于任意 n+1 个数据点,拉格朗日插值多项式可以精确地通过这些点。
  • 计算复杂性:随着数据点的增加,插值多项式的度数增加,计算复杂性呈指数增长。
  • Runge现象:在数据点较多或分布不均匀时,插值多项式可能出现振荡现象,导致插值效果变差。
比较
  • 计算复杂性:线性插值 < 三次样条插值 < 拉格朗日插值。
  • 平滑性:线性插值 < 拉格朗日插值 < 三次样条插值。
  • 适用场景
    • 线性插值适用于简单且数据变化均匀的情况。
    • 三次样条插值适用于需要平滑曲线的情况。
    • 拉格朗日插值适用于数据点较少且需要精确插值的情况。

参考资料:

Datawhale 数学建模教程:第六章

  • 27
    点赞
  • 18
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值