R语言在分析化学教学中的应用
Application of Programming Language R in Teaching of Analytical Chemistry
通讯作者:
收稿日期: 2020-05-5 接受日期: 2020-05-8
Received: 2020-05-5 Accepted: 2020-05-8
Application of programming language R for the teaching of analytical chemistry was demonstrated with several examples. Numerical calculation of R was used to solve the high-order equations involved in chemical equilibrium. It avoided the cumbersome processes of formula simplification, and also unified all equilibrium calculations into the same method. Plotting functions of R were used to plot species distribution curves, titration curves, etc. Statistical analysis functions of R were used to replace the traditional manual method. Numerical simulation functions of R were used to help students understand complex concepts such as propagation of error and confidence interval. Compared with Excel, which is commonly used in the teaching of analytical chemistry, R also has the characteristics of being code-based, self-explanatory, flexible, rapidly evolving, and readily accessible. The application of R in teaching of analytical chemistry not only helps to improve the learning efficacy, but also provides an opportunity for students to master the powerful research tool R.
Keywords:
本文引用格式
谭巧国.
Tan Qiaoguo.
1 R语言简介
R语言是一款开源免费的编程语言,主要用于统计分析和数据可视化。R最初是由两位统计学家罗斯∙伊哈卡和罗伯特∙杰特曼基于S语言开发而成。由于两位开发者名字的首字母都是R,R因此得名。目前R由“R开发核心团队”负责维护和开发[9]。
R的功能可以通过程序包得到扩展。例如本文中绘图用到了ggplot2程序包,数据整理用到了reshape2程序包,而这两个程序包与其他多个常用的程序包一起又包括在tidyverse程序包中。初次使用时,需运行以下代码安装:
install.packages("tidyverse")
除了以上提及的统计分析和绘图功能,R用于分析化学教学还具有以下几方面的优势:
R以代码形式实现计算和作图,操作过程完整地体现在代码中。代码具有自释性,理论上无需另外的解释性文字。代码可以纯文本的形式储存和传播,操作过程能够通过代码精确重现。与此相比,Excel计算过程的展示则必须用到图片(甚至视频)结合繁琐的步骤说明[11],不仅占用大量储存空间,而且也不利于教学中的传播、分享和交流。
R有广阔的应用场景。R作为以数据分析和可视化见长的编程语言,是科学研究的有力工具。R借助扩展程序包,不断发展出新的功能满足科研需求。在分析化学的教学中使用R,为学习者接触掌握这项科研利器提供良好契机,可提高学习者的知识获得感和实际收益。
R软件免费,轻便。在写作本文时,最新版本R 4.0.0的安装文件仅84 MB。因此便于在不同的场合(如教室、师生的个人电脑等)随时安装使用R。与此相比,Matlab等收费软件不仅价格昂贵,而且需要5–8 GB的安装空间。
2 R应用于分析化学教学举例
2.1 R应用于分析化学计算
在学习化学平衡和滴定分析时,常涉及复杂的计算。例如,一元弱酸溶液pH值的计算式是一元三次方程,多元弱酸和其他溶液还会涉及更高次的方程。高次方程不一定存在公式解,即使存在其公式也极为复杂。因此,直接解高次方程缺乏可操作性。目前教材中常采用忽略次要组分,忽略方程中的次要项等做法,以损失一定的准确度为代价,将高次方程简化成二次或一次方程进行求解[12]。简化中涉及多个判别式,判别式因溶液体系而异。这部分内容因此记忆负担较大,常成为学习分析化学的障碍。
R语言可看成一个功能全面的计算器,能执行各类科学计算,能用于解决分析化学中的高次方程求解问题。下面举一例说明求解过程。例题取自柳青和王海水[13]的论文,关于磷酸溶液pH值的计算。柳青和王海水提出了简化多元弱酸溶液pH值计算的新方法,但思路仍然是将高次方程简化为二次或一次方程。其方法改善了通用性,但过程仍然较为复杂,仍有较大的记忆负担。
例1 已知磷酸(H3A)的Ka1 = 7.6 × 10−3,Ka2 = 6.3 × 10−8,Ka3= 4.4 × 10−13,计算0.20 mol∙L−1磷酸溶液的pH。
以上例题可通过执行以下R代码进行计算。需说明的是#号后是代码注释,用于说明代码的用途,注释不会被作为代码执行。本文中的代码在关键处添加了注释,便于读者理解。
#输入已知条件
Ka1 = 7.6E-3
Ka2 = 6.3E-8
Ka3 = 4.4E-13
c = 0.20
#定义三元酸H+浓度的精确计算式
f = function(H) {
#将分母命名为D
D = H^3 + H^2*Ka1 + H*Ka1*Ka2 + Ka1*Ka2*Ka3
H - 1E-14/H - (c*H^2*Ka1 + c*2*H*Ka1*Ka2 + c*3*Ka1*Ka2*Ka3)/D
}
#运用uniroot函数解方程
H = uniroot(f, interval = c(0, 0.6), tol = 10^-20)$root
#根据H+浓度计算
pH = -log10(H)
pH
以上代码中,起核心作用的是解方程函数uniroot()。在使用uniroot()时,所输入的f是H+浓度精确计算式,因溶液体系而异;参数interval是粗略预估的方程根的上下限;参数tol是所允许的误差。运行以上代码后,得到结果:1.45134,因此可知溶液pH值为1.45,与参考文献中的计算值一致。
运用该方法,除了计算准确,还可以将所有酸碱溶液的pH值计算都转变为同一种方法。计算过程可分为三步:(1)根据质子得失平衡写出质子条件式;(2)将式中酸(碱)各型体的浓度都替换成酸(碱)的分析浓度与分布分数之积,得到关于[H+]([OH−])的一元高次方程,即以上所用的f;(3)解高次方程。
此方法只用到质子条件式和分布分数这两个容易掌握的基本方法,避免了复杂的判别式和简化过程,便于理解和记忆。从而使学习者从复杂的公式和繁琐的计算中解脱出来,将更多精力用于学习分析化学知识本身。
2.2 R应用于分析化学绘图
R语言另一项特色功能是科学绘图。图在分析化学的教学中发挥着不可或缺的作用,可以将复杂的概念、溶液体系、反应过程等可视化,降低理解难度。例如,分布分数图之于化学平衡,滴定曲线图之于滴定原理,都是极为有效的教学工具。学习者若能自己绘制这些图,则更能提升学习效果,增加学习乐趣。
以下以多元酸的分布分数图和络合滴定图为例,演示如何用将R应用于分析化学绘图,并帮助理解其背后的化学原理。需要说明的是,为了节省篇幅,以下作图的R代码仅保留主干部分,图的细节修饰所需代码未给出。读者运行以下代码,可看到图的主体信息,但格式细节与文中所示的图略有区别。完整作图代码获取方式见文末说明。
2.2.1 分布分数图的绘制
三元弱酸磷酸经三级解离,得到四种型体。根据磷酸的三个Ka值,可得到四种型体分布分数的表达式,并据此绘制分布分数图。以下是绘图的R代码:
library(ggplot2) #加载用于绘图的程序包
Ka1 = 7.6E-3
Ka2 = 6.3E-8
Ka3 = 4.4E-13
#选取0~14范围内的100个pH值
pH = seq(from = 0, to = 14, length.out = 100)
H = 10^(-pH)
#根据H+浓度计算分布分数的分母
D = H^3 + H^2*Ka1 + H*Ka1*Ka2 + Ka1*Ka2*Ka3
#计算磷酸各型体的分布分数
delta_0 = Ka1*Ka2*Ka3/D
delta_1 = H*Ka1*Ka2/D
delta_2 = H^2*Ka1/D
delta_3 = H^3/D
#将pH值和各型体分布分数组合成数据表,用于绘图
df = data.frame(pH, delta_0, delta_1, delta_2, delta_3)
#绘制分布分数图
ggplot(df, aes(x = pH))+ #x轴为pH
geom_line(aes(y = delta_0), linetype = 1)+ #PO4对应的曲线
geom_line(aes(y = delta_1), linetype = 2)+ #HPO4对应的曲线
geom_line(aes(y = delta_2), linetype = 3)+ #H2PO4对应的曲线
geom_line(aes(y = delta_3), linetype = 4) #H3PO4对应的曲线
运行以上代码,可得到磷酸的分布分数图(图 1)。该图可用于判断不同pH值条件下,磷酸的主要型体。例如,天然海水的pH约为8.2,因此可判断海水中的磷酸盐主要以
图1
除了一元、多元酸碱的分布分数图,其他分布分数图,如Cu(NH3)n逐级络合物的分布分数图,EDTA各型体的分布分数图等,亦可按照以上方法绘制。学习者若改变作图参数(如Ka值),便可即时看到图的变化。这样的过程可以帮助理解稳定常数、解离常数等如何决定了曲线的形状,进而理解这些常数的内涵。
2.2.2 滴定曲线的绘制
以下以EDTA络合滴定金属离子M为例,演示如何用R绘制络合滴定曲线。络合滴定曲线是滴定过程中金属离子浓度(pM')与滴定分数(a)之间关系的曲线,由络合滴定曲线方程描述:
在此例中,假设用0.1000 mol∙L−1的EDTA滴定同浓度的金属离子M,EDTA络合M的条件稳定常数K'MY为1.0 × 1014 L∙mol−1。绘制其滴定曲线图的R代码如下:
#设定金属离子和EDTA溶液的分析浓度
cM = 0.1
#设定条件稳定常数值
KMY = 1E14
#设定滴定分数,在0~2之间取200个点
a_seq = seq(0, 2, length.out = 200)
#设定名为pM,长度为200的向量,用于储存计算所得的pM值
pM = numeric(200)
#针对每一个滴定分数值a,解络合滴定方程计算对应的pM值,此处运用循环语句
for (i in 1:200) {
#定义络合滴定曲线方程
f = function(x) {
#将已知条件导入方程
a < < - a_seq[i]
KMY < < - KMY
cM < < - cM
KMY*x^2+(KMY*cM/(1+a)*(a-1)+1)*x-cM/(1+a)
}
#将方程的根存入名为M的变量
M = uniroot(f, interval = c(1E-20, 0.1), tol=1E-20)$root
#计算pM值并保存
pM[i] = -log10(M)
}
#将滴定分数值与对应的pM值合并成数据表,命名为df
df = data.frame(a = a_seq, pM = pM)
#调用作图所需软件包ggplot2
library(ggplot2)
#作图
ggplot(df, aes(x = a, y = pM))+
geom_line()+
labs(x = "滴定分数", y = "pM'")
运行以上代码,得到络合滴定曲线(图 2)。曲线中包含了络合滴定的多项信息,这些信息在曲线绘制过程中均得到求解。例如,化学计量点的pM'即滴定分数a等于1时的pM',为7.65;滴定突跃范围是滴定分数a为(100.0 ± 0.1)%对应的pM'范围,为4.30–11.00。通过调整K'MY和cM的取值,观察曲线形状的变化,有助于理解络合稳定常数以及金属离子浓度对于滴定突跃的影响。
图2
2.3 R应用于分析化学中的统计分析
定量分析中的误差不可避免,计算和评估误差则需要用到统计分析方法。统计分析的理论和方法较为抽象,对其理解有赖于扎实的数学和统计学基础,因此常成为教和学的难点。分析化学中统计分析的教学应有别于数学、统计学类课程,目标应当是准确理解,正确使用,可以淡化数学推导。R最初是由统计学家开发,其统计学功能齐全而严谨,且容易使用。将R运用于分析化学中统计方法的教学,不仅可以帮助学习者掌握数据处理的正确方法,还有助于更好地理解相关的理论和方法。
以下以t检验为例,演示如何将R用于分析化学中的统计分析,并结合可视化的方法帮助学习者理解统计学概念。例题取自武汉大学《分析化学(上册)》[12]。
例2 采用一种新方法测定基准明矾中铝的质量分数,9次测定结果为10.74%,10.77%,10.77%,10.77%,10.81%,10.82%,10.73%,10.86%,10.81%。已知明矾中铝含量的标准值(以理论值代替)为10.77%。采用新方法后,是否引起系统误差(置信度95%)?
教材上介绍的是常规方法:由于测定次数较少,判断随机误差符合t分布。计算t值得到1.43,查t值表得到临界值t(0.05, 8) = 2.31。t值小于临界值,从而判断无显著系统误差。
运用R语言,无需查表,进行单样本t检验直接得到结果。具体操作如下:
#将9次测定值存入变量x
x = c(10.74, 10.77, 10.77, 10.77, 10.81, 10.82, 10.73, 10.86, 10.81)
#用t.test函数进行t检验,参数mu是理论值,参数conf.level是置信度
t.test(x, mu = 10.77, conf.level = 0.95)
输出结果如下(为节省篇幅,部分信息省略):
## One Sample t-test
## t = 1.2039, df = 8, p-value = 0.2631
## alternative hypothesis: true mean is not equal to 10.77
## 95 percent confidence interval:
## 10.75474 10.81859
结果显示p值为0.2631,大于显著性水平0.05,因而可判断新方法的测定值与理论值无显著差异,亦即新方法未导致显著的系统误差。结果中还给出了测定值的95%置信区间,即(10.755,10.819),此区间内包含理论值10.77,也表明新方法无显著系统误差。
图3
图3
用R绘图对比测定值及其均值的95%置信区间与理论值
(a) 95%置信区间包含理论值,无显著系统误差;(b) 95%置信区间不包含理论值,存在显著系统误差;(c) 99%置信区间比95%置信区间宽,包含理论值,无显著系统误差
library(ggplot2)
ggplot(data.frame(x), aes(x, y=0))+
geom_rect(xmin=10.75474, xmax=10.81859, ymin=-0.01, ymax=0.01)+ #置信区间
geom_point(alpha = 1/3)+ #各测定值对应的点
geom_text(label=x, nudge_y=0.01)+ #各测定值
geom_point(x=10.77, y=0, shape= 1, size =4)+ #理论值
ylim(-0.1, 0.1)
从图上可看到各次测定值的位置决定了置信区间的位置和宽度;置信区间若将理论值包含在内,则表明测定值的均值与理论值无显著差异。为了进一步理解,另外生成一组有系统误差的数据作为对比。将例题中的测定值均减去0.05,并针对这组新数据做t检验:
x_2 = x - 0.05 #原数据各减去0.05
t.test(x_2, mu = 10.77, conf.level = 0.95)
结果显示p值为0.04266,小于0.05;95%置信区间为(10.705, 10.768),不包含理论值10.77。表明存在显著的系统误差(图 3b)。若将显著性水平调为99%,再进行t检验:
t.test(x_2, mu = 10.77, conf.level = 0.99)
得到结果p值仍为0.04266,虽小于0.05,但大于新采用的显著性水平0.01;99%置信区间为(10.690, 10.783),比95%置信区间更宽,包含理论值10.77,表明不存在显著系统误差。因此可知,关于是否存在系统误差的结论与置信度的选择有密切相关。
除以上的t检验外,分析化学教学中涉及正态分布、t分布、置信区间的估计、F检验、可疑值的取舍、回归分析等统计分析内容均可借助R进行教学。
2.4 R数值模拟用于解释复杂概念
R作为编程语言,有生成符合特定分布的随机数,进行随机抽样等功能,因而可以进行数值模拟,帮助解释分析化学中一些复杂的概念和方法,以下举一例说明。
例3 测量值m1 = 0.00 ± 0.10,测量值m2 = 2.00 ± 0.10,则m2与m1之差m的值和标准差是多少?
按照误差传递公式,可以算得m的标准差是(0.102 + 0.102)1/2,即0.14,因此m = 2.00 ± 0.14。误差传递计算并不难,但是误差传递公式的推导过程较为复杂,学习和理解的难度较高。因此,学习者虽然能够按照公式算出误差,但未必能理解误差传递的真正含义。借助R的数值模拟和可视化,可以帮助学习者克服这一困难。以下针对例题3,进行数值模拟,得到图 4。
图4
library(ggplot2)
library(reshape2)
#随机生成4组包含100个符合正态分布的数,参数mean为均值,sd为标准偏差
m1 = rnorm(100, mean = 0, sd = 0.10) #均值为0,标准差为0.10
m2 = rnorm(100, mean = 2, sd = 0.10) #均值为2,标准差为0.10
r1 = rnorm(100, mean = 2, sd = 0.20) #均值为2,标准差为0.20
r2 = rnorm(100, mean = 2, sd = 0.1414) #均值为2,标准差为0.1414
#计算m2与m1之差
m = m2 - m1
#运用data.frame函数将上述5组数合组合成数据表d
d = data.frame(m1, m2, r1, r2, m)
#运用gather函数排列的数据表d的转化成纵向排列的数据表d_1
d1 < - melt(d)
#运用ggplot函数作图
ggplot(d1, aes(x = variable, y = value))+
geom_violin()+
geom_jitter(width = 0.1, shape = 1)
3 结语
以上举例展示了R在分析化学教学中的应用。R兼具计算、绘图、统计、模拟等功能,因此,分析化学中凡是涉及这些功能的内容,不限于以上所举的例子,均可用R辅助教学。R通过代码实现各项功能,操作过程均体现在文本化的代码中,便于信息的传播交流。因此R不仅适用于课堂演示,还便于学习者课后自学。相比于Excel,R在初学时有较为陡峭的学习曲线。考虑到本科生通常在低年级学习过C语言等更为复杂的编程语言,学习R已经有了较好的基础。目前在高校和学术界有人数众多气氛活跃的R学习群体,互联网上有丰富且优质的R学习资料,这也为自学R营造了良好的环境。将R用于分析化学的教与学,也让学习者有机会在应用实践中学习R,反过来促进R的学习。从分析化学课上学到R的技能,也将有助于其他课程的学习和日后的科学研究。
本文中计算、绘图等所涉及R代码的完整版本均可在此链接(http://tinyurl.com/y9hh36zq)下载或通过联系作者获取。
参考文献
/
〈 |
|
〉 |
