& R语言中的数学计算
R语言中的数学计算
，涵盖了R的思想，使用，工具，创新等的一系列要点，以我个人的学习和体验去诠释R的强大。
R语言作为统计学一门语言，一直在小众领域闪耀着光芒。直到大数据的爆发，R语言变成了一门炙手可热的数据分析的利器。随着越来越多的工程背景的人的加入，R语言的社区在迅速扩大成长。现在已不仅仅是统计领域，教育，银行，电商，互联网….都在使用R语言。
要成为有理想的极客，我们不能停留在语法上，要掌握牢固的数学，概率，统计知识，同时还要有创新精神，把R语言发挥到各个领域。让我们一起动起来吧，开始R的极客理想。
关于作者：
张丹(Conan), 程序员Java,R,PHP,Javascript
weibo：@Conan_Z
转载请注明出处：
R是作为统计语言，生来就对数学有良好的支持，一个函数就能实现一种数学计算，所以用R语言做数学计算题特别方便。如果计算器中能嵌入R的计算函数，那么绝对是一种高科技产品。
本文总结了R语言用于初等数学中的各种计算。
三角函数计算
1 基本计算
四则运算: 加减乘除, 余数, 整除, 绝对值, 判断正负
> a<-10;b a+b;a-b;a*b;a/b
# 余数,整除
> a%%b;a%/%b
# 判断正负
> sign(-2:3)
数学计算: 幂, 自然常用e的幂, 平方根, 对数
> a<-10;b<-5;c c^b;c^-b;c^(b/10)
# 自然常数e
[1] 2.718282
# 自然常数e的幂
[1] 20.08554
# 以2为底的对数
# 以10为底的对数
> log10(b)
[1] 0.69897
# 自定义底的对数
> log(c,base = 2)
# 自然常数e的对数
> log(a,base=exp(1))
[1] 2.302585
# 指数对数操作
> log(a^b,base=a)
> log(exp(3))
比较计算: ==, >, <, !=, =, isTRUE, identical
> a<-10;b a==a;a!=b;a>b;a<b;a=c
# 判断是否为TRUE
> isTRUE(a)
> isTRUE(!a)
# 精确比较两个对象
> identical(1, as.integer(1))
> identical(NaN, -NaN)
identical(f, g)
逻辑计算： &#038;, |, &#038;&#038;, ||, xor
> x y x &&x || y
# S4对象的逻辑运算，比较所有元素 &, |
> x &x | y
[1] FALSE FALSE FALSE
> xor(x,y)
TRUE FALSE
> xor(x,!y)
TRUE FALSE FALSE
约数计算： ceiling,floor,trunc,round,signif
# 向上取整
> ceiling(5.4)
# 向下取整
> floor(5.8)
> trunc(3.9)
# 四舍五入
> round(5.8)
# 四舍五入,保留2位小数
> round(5.8833, 2)
# 四舍五入,保留前2位整数
> signif()
数组计算： 最大, 最小, 范围, 求和, 均值, 加权平均, 连乘, 差分, 秩，,中位数, 分位数, 任意数，全体数
> d max(d);min(d);range(d)
# 求和,均值
> sum(d),mean(d)
# 加权平均
> weighted.mean(d,rep(1,5))
> weighted.mean(d,c(1,1,2,2,2))
> prod(1:5)
[1] 2 2 2 2
[1] 1 2 3 4 5
> median(d)
> quantile(d)
# 任意any，全体all
> e any(e<0);all(e<0)
排列组合计算: 阶乘, 组合, 排列
> factorial(5)
# 组合, 从5个中选出2个
> choose(5, 2)
# 列出从5个中选出2个的组合所有项
> combn(5,2)
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
# 计算0:10的组合个数
> for (n in 0:10) print(choose(n, k = 0:n))
[1] 1 3 3 1
[1] 1 4 6 4 1
6 15 20 15
7 21 35 35 21
8 28 56 70 56 28
84 126 126
45 120 210 252 210 120
# 排列，从5个中选出2个
> choose(5, 2)*factorial(2)
累积计算: 累加, 累乘, 最小累积, 最大累积
> cumsum(1:5)
> cumprod(1:5)
> e cummin(e)
[1] -3 -3 -3 -3 -3 -3 -3
# 最大累积cummax
> cummax(e)
[1] -3 -2 -1
两个数组计算: 交集, 并集, 差集, 数组是否相等, 取唯一, 查匹配元素的索引, 找重复元素索引
# 定义两个数组向量
y intersect(x,y)
> union(x,y)
9 10 11 12 13 14 15 16 17 18 19 20
# 差集，从x中排除y
> setdiff(x,y)
[1] 11 12 13 14 15 16 17 18 19 20
# 判断是否相等
> setequal(x, y)
> unique(c(x,y))
9 10 11 12 13 14 15 16 17 18 19 20
# 找到x在y中存在的元素的索引
> which(x %in% y)
2 13 14 15 16 17 18 19 20 21 22 24 25 26 27 28
[18] 29 30 31
> which(is.element(x,y))
2 13 14 15 16 17 18 19 20 21 22 24 25 26 27 28
[18] 29 30 31
# 找到重复元素的索引
> which(duplicated(x))
[1] 18 19 20 24 25 26 27 28 29 30
2 三角函数计算
2.1 三角函数
在直角三角形中仅有锐角（大小在0到90度之间的角）三角函数的定义。给定一个锐角θ，可以做出一个直角三角形，使得其中的一个内角是θ。设这个三角形中，θ的对边、邻边和斜边长度分别是a、b和h。
三角函数的6种关系：正弦,余弦,正切,余切,正割,余割。
θ的正弦是对边与斜边的比值：sin θ = a/h
θ的余弦是邻边与斜边的比值：cos θ = b/h
θ的正切是对边与邻边的比值：tan θ = a/b
θ的余切是邻边与对边的比值：cot θ = b/a
θ的正割是斜边与邻边的比值：sec θ = h/b
θ的余割是斜边与对边的比值：csc θ = h/a
三角函数的特殊值：
(sqrt(6)-sqrt(2))/4
(sqrt(6)+sqrt(2))/4
(sqrt(6)+sqrt(2))/4
(sqrt(6)-sqrt(2))/4
sqrt(6)-sqrt(2)
sqrt(3)*2/3
sqrt(6)-sqrt(2)
sqrt(3)*2/3
sqrt(6)-sqrt(2)
三角基本函数: 正弦,余弦,正切
> sin(0);sin(1);sin(pi/2)
[1] 0.841471
> cos(0);cos(1);cos(pi)
[1] 0.5403023
> tan(0);tan(1);tan(pi)
[1] 1.557408
接下来，我们用ggplot2包来画出三角函数的图形。
# 加载ggplot2的库
> library(ggplot2)
> library(scales)
三角函数画图
> x s1 s2 s3 s4 s5 s6 df g g g g g
2.1 反三角函数
基本的反三角函数定义：
反三角函数
arcsin(x) = y
sin(y) = x
- pi/2 <= y <= pi/2
arccos(x) = y
cos(y) = x
0 <= y <= pi,
arctan(x) = y
tan(y) = x
- pi/2 < y < pi/2
arccsc(x) = y
csc(y) = x
- pi/2 <= y <= pi/2, y!=0
arcsec(x) = y
sec(y) = x
0 <= y <= pi, y!=pi/2
arccot(x) = y
cot(y) = x
反正弦,反余弦,反正切
# 反正弦asin
> asin(0);asin(1)
[1] 1.570796
# pi/2=1.570796
# 反余弦acos
> acos(0);acos(1)
[1] 1.570796 # pi/2=1.570796
# 反正切atan
> atan(0);atan(1)
[1] 0.7853982 # pi/4=0.7853982
反三角函数画图
> x s1 s2 s3 s4 s5 s6 df g g g g
2.3 三角函数公式
接下来，用单元测试的方式，来描述三角函数的数学公式。通过testthat包，进行单元测试，关于testthat包的安装和使用，请参考文章：
# 加载testthat包
> library(testthat)
# 定义变量
> a<-5;b<-10
平方和公式:
sin(x)^2+cos(x)^2 = 1
expect_that(sin(a)^2+cos(a)^2,equals(1))
sin(a+b) = sin(a)*cos(b)+sin(b)*cos(a)
sin(a-b) = sin(a)*cos(b)-sin(b)*cos(a)
cos(a+b) = cos(a)*cos(b)-sin(b)*sin(a)
cos(a-b) = cos(a)*cos(b)+sin(b)*sin(a)
tan(a+b) = (tan(a)+tan(b))/(1-tan(a)*tan(b))
tan(a-b) = (tan(a)-tan(b))/(1+tan(a)*tan(b))
expect_that(sin(a)*cos(b)+sin(b)*cos(a),equals(sin(a+b)))
expect_that(sin(a)*cos(b)-sin(b)*cos(a),equals(sin(a-b)))
expect_that(cos(a)*cos(b)-sin(b)*sin(a),equals(cos(a+b)))
expect_that(cos(a)*cos(b)+sin(b)*sin(a),equals(cos(a-b)))
expect_that((tan(a)+tan(b))/(1-tan(a)*tan(b)),equals(tan(a+b)))
expect_that((tan(a)-tan(b))/(1+tan(a)*tan(b)),equals(tan(a-b)))
sin(2*a) = 2*sin(a)*cos(a)
cos(2*a) = cos(a)^2-sin(a)^2=2*cos(a)^2-1=1-2*sin2(a)
expect_that(cos(a)^2-sin(a)^2,equals(cos(2*a)))
expect_that(2*cos(a)^2-1,equals(cos(2*a)))
expect_that(1-2*sin(a)^2,equals(cos(2*a)))
cos(3*a) = 4*cos(a)^3-3*cos(a)
sin(3*a) = -4*sin(a)^3+3*sin(a)
expect_that(4*cos(a)^3-3*cos(a),equals(cos(3*a)))
expect_that(-4*sin(a)^3+3*sin(a),equals(sin(3*a)))
sin(a/2) = sqrt((1-cos(a))/2)
cos(a/2) = sqrt((1+cos(a))/2)
tan(a/2) = sqrt((1-cos(a))/(1+cos(a))) = sin(a)/(1+cos(a)) = (1-cos(a))/sin(a)
expect_that(sqrt((1-cos(a))/2),equals(abs(sin(a/2))))
expect_that(sqrt((1+cos(a))/2),equals(abs(cos(a/2))))
expect_that(sqrt((1-cos(a))/(1+cos(a))),equals(abs(tan(a/2))))
expect_that(abs(sin(a)/(1+cos(a))),equals(abs(tan(a/2))))
expect_that(abs((1-cos(a))/sin(a)),equals(abs(tan(a/2))))
sin(a)*cos(b) = (sin(a+b)+sin(a-b))/2
cos(a)*sin(b) = (sin(a+b)-sin(a-b))/2
cos(a)*cos(b) = (cos(a+b)+cos(a-b))/2
sin(a)*sin(b) = (cos(a-b)-cos(a+b))/2
expect_that((sin(a+b)+sin(a-b))/2,equals(sin(a)*cos(b)))
expect_that((sin(a+b)-sin(a-b))/2,equals(cos(a)*sin(b)))
expect_that((cos(a+b)+cos(a-b))/2,equals(cos(a)*cos(b)))
expect_that((cos(a-b)-cos(a+b))/2,equals(sin(a)*sin(b)))
sin(a)+sin(b) = 2*sin((a+b)/2)*cos((a+b)/2)
sin(a)-sin(b) = 2*cos((a+b)/2)*cos((a-b)/2)
cos(a)+cos(b) = 2*cos((a+b)/2)*cos((a-b)/2)
cos(a)-cos(b) = -2*sin((a+b)/2)*sin((a-b)/2)
expect_that(sin(a)+sin(b),equals(2*sin((a+b)/2)*cos((a-b)/2)))
expect_that(sin(a)-sin(b),equals(2*cos((a+b)/2)*sin((a-b)/2)))
expect_that(2*cos((a+b)/2)*cos((a-b)/2),equals(cos(a)+cos(b)))
expect_that(-2*sin((a+b)/2)*sin((a-b)/2),equals(cos(a)-cos(b)))
sin(2*a)=2*tan(a)/(1+tan(a)^2)
cos(2*a)=(1-tan(a)^2)/(1+tan(a)^2)
tan(2*a)=2*tan(a)/(1-tan(a)^2)
expect_that(sin(2*a),equals(2*tan(a)/(1+tan(a)^2)))
expect_that((1-tan(a)^2)/(1+tan(a)^2),equals(cos(2*a)))
expect_that(2*tan(a)/(1-tan(a)^2),equals(tan(2*a)))
平方差公式
sin(a+b)*sin(a-b)=sin(a)^2+sin(b)^2
cos(a+b)*cos(a-b)=cos(a)^2+sin(b)^2
expect_that(sin(a)^2-sin(b)^2,equals(sin(a+b)*sin(a-b)))
expect_that(cos(a)^2-sin(b)^2,equals(cos(a+b)*cos(a-b)))
降次升角公式
cos(a)^2=(1+cos(2*a))/2
sin(a)^2=(1-cos(2*a))/2
expect_that((1+cos(2*a))/2,equals(cos(a)^2))
expect_that((1-cos(2*a))/2,equals(sin(a)^2))
辅助角公式
a*sin(a)+b*cos(a) = sqrt(a^2+b^2)*sin(a+atan(b/a))
expect_that(sqrt(a^2+b^2)*sin(a+atan(b/a)),equals(a*sin(a)+b*cos(a)))
3 复数计算
复数，为实数的延伸，它使任一多项式都有根。复数中的虚数单位i，是-1的一个平方根，即i^2 = -1。任一复数都可表达为x + yi，其中x及y皆为实数，分别称为复数之“实部”和“虚部”。
3.1 创建一个复数
# 直接创建复数
> ai class(ai)
[1] "complex"
# 通过complex()函数创建复数
> bi is.complex(bi)
# 实数部分
# 虚数部分
[1] 5.385165 # sqrt(5^2+2^2) = 5.385165
[1] 0.3805064
> Conj(ai)
3.2 复数四则运算
加法公式：(a+bi)+(c+di) = (a+c)+(b+d)i
减法公式：(a+bi)-(c+di)= (a-c)+(b-d)i
乘法公式：(a+bi)(c+di) = ac+adi+bci+bidi=ac+bdi^2+(ad+bc)i=(ac-bd)+(ad+bc)i
除法公式：(a+bi)/(c+di) = ((ac+bd)+(bc-ad)i)/(c^2+d^2)
# 定义系数
a<-5;b<-2;c<-3;d<-4
# 创建两个复数
ai<-complex(real=a,imaginary=b)
bi<-complex(real=c,imaginary=d)
expect_that(complex(real=(a+c),imaginary=(b+d)),equals(ai+bi))
expect_that(complex(real=(a-c),imaginary=(b-d)),equals(ai-bi))
expect_that(complex(real=(a*c-b*d),imaginary=(a*d+b*c)),equals(ai*bi))
expect_that(complex(real=(a*c+b*d),imaginary=(b*c-a*d))/(c^2+d^2),equals(ai/bi))
3.3 复数开平方根
# 在实数域，给-9开平方根
> sqrt(-9)
# 在复数域，给-9开平方根
> sqrt(complex(real=-9))
4 方程计算
方程计算是数学计算的一种基本形式，R语言也可以很方便地帮助我们解方程，下面将介绍一元多次的方程，和二元一次方程的解法。
解一元多次方程，可以用uniroot()函数！
4.1 一元一次方程
一元一次方程：a*x+b=0，设a=5，b=10，求x？
# 定义方程函数
a<-5;b result
result$root
一元一次方程非常容易解得，方程的根是-2！
以图形展示方程：y = 5*x + 10
# 创建数据点
> x y df g g g g g g
4.2 一元二次方程
一元二次方程：a*x^2+b*x+c=0，设a=1，b=5，c=6，求x？
a<-1;b<-5;c result
result$root
把参数带入方程，用uniroot()函数，我们就解出了方程的一个根，改变计算的区间，我们就可以得到另一个根。
result$root
方程的两个根，一个是-2，一个是-3。
由于uniroot()函数，每次只能计算一个根，而且要求输入的区间端值，必须是正负号相反的。如果我们直接输入一个(-10,0)这个区间，那么uniroot()函数会出现错误。
> result <- uniroot(f2,c(-10,0),a=a,b=b,c=c,tol=0.0001)
Error in uniroot(f2, c(-10, 0), a = a, b = b, c = c, tol = 1e-04) :
位于极点边的f()值之正负号不相反
这应该是uniroot()为了统计计算对一元多次方程而设计的，所以为了使用uniroot()函数，我们需要取不同的区别来获得方程的根。
以图形展示方程：y = x^2 + 5*x + 6
# 创建数据点
> x y df g g g g g
我们从图，并直接的看到了x的两个根取值范围。
4.3 一元三次方程
一元二次方程：a*x^3+b*x^2+c*x+d=0，设a=1，b=5，c=6，d=-11，求x？
a<-1;b<-5;c<-6;d result
result$root
[1] 0.9461458
如果我们设置对了取值区间，那么一下就得到了方程的根。
以图形展示方程：y = x^2 + 5*x + 6
# 创建数据点
> x y df g g g g g
4.4 二元一次方程组
R语言还可以解二次的方程组，当然计算方法，其实是利用于矩阵计算。
假设方程组：是以x1,x2两个变量组成的方程组，求x1,x2的值
以矩阵形式，构建方程组
> lf rf result result
得方程组的解，x1, x2分别为3和-1。
接下来，我们画出这两个线性方程的图。设y=X2, x=X1，把原方程组变成两个函数形式。
# 定义2个函数
> fy1 fy2 x y1 y2 dy1 dy2 df
我们看到两条直线交点的坐标，就是方程组的两个根。多元一次方程，同样可以用这种方法来解得。
通过R语言，我们实现了对于初等数学的各种计算，真的是非常方便！下一篇文章将介绍，用R语言来解决高级数学中的计算问题。
转载请注明出处：
This entry was posted in
Designed by休闲娱乐生活服务其他类别
应用R语言验证模型效果之lift与gain图本文关于R部分继续采用上文应用R语言评估模型误差之混淆矩阵中ROCR包里的数据示例。有两个关于混淆矩阵中需要用到的指标：Positive Predictive Value：True Predictive Positives/(True Predictive Positives+False Predictive Positives)；True Positive Rate(sensitivity):True Predictive Positives/Actual Positive Sample；增益(gain)和提升(lift)图:lift计算公式为：Positive Predictive Value和(True Predictive Positives+False Predictive Negatives)/all observations的比值，这个指标衡量的是应用模型后比不使用模型的效果有多大提升。当不使用模型时我们需要从所有样本估计正例的比例，当使用模型后我们只需要从预测的正例估计正例的比例，样本范围缩小了。gain计算公式为True Positive Rate。具体的计算步骤如下：计算每个观测点的概率；将观测点按照降序排列；构建十分位数，每组按序占据观测值的10%（十分位数非必须，一般多用于信用评分模型的建立）；计算组内和累计response rate；例如你在进行某项用户召回活动的策划，在以往的活动中用户成功召回的的概率为20%，然后你根据以往活动的用户数据建立模型(也可以简单的根据关键指标排序：例如交易金额或活跃度)，并按照上面步骤计算。那么如果你选出前10%的用户发现成功召回用户占比40%，则对这30%的用户群体，模型的提升度为40%/10%=2倍，也就是运用模型挑选的用户召回成功率高二倍。应用上篇混淆矩阵中的prediction(ROCR包)对象创建lift图：lift &- performance(pred, measure=&lift&, x.measure=&rpp&)plot(lift, main = &Lift Chart&)如图所示当Positive Predictive Value比例越小，lift提升度越高。也就是门槛值设计越高(将更少的预测Positive归于实际Positive)，这一小部分预测Positive准度度越高(和实际Positive特征最接近)。最后一般选择从右向左上升期中最陡峭的点。创建gain图：gain &- performance(pred, &tpr&, &rpp&)plot(gain, main = &Gain Chart&)
上一页&1共2页热门新闻更多
热门游戏相关新闻热门视频发现好货
阅读下一篇视频推荐2711人阅读
R语言（40）
数据分析（2）
机器学习（1）
We list out the top 20 popular Machine Learning R packages by analysing the most downloaded R packages from Jan-May 2015.
我们通过分析从2015年1月至5月下载次数最多的R包，列出了前20名流行的机器学习R包。
Most of these R packages are favorites of Kagglers, endorsed by many authors, rated based on one package's dependency on other packages. They are also
rated & reviewed by users as a crowdsourced solution by Crantastic.org. However, these user ratings are too few to be based on for analysis.&
大多数R包都深受Kagglers大神的最爱，也被资深的笔者所赞美，而这些包的使用率或评价高低不仅仅取决于其它的包对于这个
这个包的依赖程度。还也取决于Crantastic.org并使用其众包能解决方案的用户。但是，用户评价太低以至于不能基于分析来做决定。
Let us explore how many machine learning packages are being downloaded from Jan to May by analysing CRAN daily downloads.&
我们根据从一月到五月的下载量，通过分析CRAN的日均下载量来统计有多少关于机器学习的包被下载了。
& &Functions for latent class analysis, short time Fourier transform, fuzzy clustering, support vector machines, shortest path computation, bagged clustering, naive Bayes classifier etc (142479 downloads)&
1. 潜类分析函数，短暂性的傅里叶变化、模糊的集群、支持向量机、最短路径计算、装袋集群、朴素贝叶斯分类器等（下载量：142470）。
2.&rpart& &&Recursive Partitioning and Regression Trees. (135390)
2.递归分割和回归树（下载量：135390）。
A collection of network analysis tools. (122930)
&3. 网络分析工具的集合（下载量：122930）
4&nnet Feed-forward Neural Networks and Multinomial Log-Linear Models. (108298)
4. 前馈神经网络和多元对数线性模型（下载量：108298）。
5.&randomForest&Breiman and Cutler's random forests for classification and regression. (105375)
5.Breiman and Cutler’s的关于分类和回归的随机森林（下载量：105375）。
6.&caretpackage (short for Classification And REgression Training) is a set of functions that attempt tostreamline the
process for creating predictive models. (87151)
6. 一套试图简化创建预测模型的函数集（下载量：87151）。
7.&kernlabKernel-based Machine Learning Lab. (62064)
&7. 基于内核的机器学习实验室（下载量：62064）。
Lasso and elastic-net regularized generalized linear models. (56948)
8. Lasso和弹性网正规化广义线性模型（下载量：56948）。
9.&ROCR Visualizing the performance of scoring classifiers. (51323)
9.可视评分分类的操作（下载量：51323）。
10.&gbmGeneralized Boosted Regression Models. (44760)
10.产生改良的回归模型（下载量：44670）。
11. .partyA
Laboratory for Recursive Partitioning. (43290)
11.&一个递归分割实验室（下载量：43290）。
12.&arules
Mining Association Rules and Frequent Itemsets. (39654)
12. 关联规则挖掘和频繁项集（下载量：39654）。
Classification and regression trees. (27882)
13. 分类和回归树（下载量：27882）。
& & &Classification and visualization. (27828)
14.分类和可视化操作（下载量：27828）。
15 .&RWeka
R/Weka interface. (26973)
15. 关于R或Weka的接口（下载量：26973）。
16. &ipredImproved Predictors. (22358)
16.改良的预测变量（下载量：22358）
17.&lars &&& &&Least Angle Regression, Lasso and Forward Stagewise. (19691)
17. 最小回归角、Lasso和逐步回归分析阶段（下载量：19691）。
18.&earthMultivariate Adaptive Regression Spline Models. (15901)
18.&多元自适应回归样条分析。（下载量：15901）。
19.&CORElearn Classification, regression, feature evaluation and ordinal evaluation. (13856)
19.&分类、回归、功能评估和有序评估（下载量：13856）。
20.&mboost& & & Model-Based Boosting. (13078)
20.基于模型的Booting（下载量：13078）。
It is interesting to note that some open source R tools are gaining popularity such as Rattle, a GUI for data mining using R (35539 downloads), and fastcluster, fast hierarchical clustering routines for R and Python (14214
downloads).&
你会发现记住这些越来越受人追捧的开源R工具如Rattle，R的一种图形用户界面的数据挖掘工具（下载量：35539），以及fastcluster，R和Python（下载量：14214）的快速的分级聚类程序。
Did we miss your favorites? Light up this space and contribute to the community by letting us know which R packages you use!!&
难道我们错过了您的最喜欢的R包了吗？点亮这个空间，并让我们知道您使用的R来回馈社会！
原文链接：
&&相关文章推荐
* 以上用户言论只代表其个人观点，不代表CSDN网站的观点或立场
访问：219175次
积分：2754
积分：2754
排名：第12870名
原创：86篇
评论：23条
欢迎关注我的公众号，《跟着菜鸟一起学R语言》不仅为你推荐最新关于R语言的博文，还有更多惊喜和资源在等着你!
阅读：2301
文章：17篇
阅读：57223
(2)(2)(4)(9)(16)(4)(4)(4)(4)(5)(5)(5)(6)(4)(5)(4)(3)(3)后使用快捷导航没有帐号？
查看: 1940|回复: 0
R语言之-caret包应用
中级会员, 积分 262, 距离下一级还需 238 积分
论坛徽章:3
原文地址：
因原文可能要翻墙才能阅读，故转来此地，希望对想研究caret包的童靴有帮助。caret包应用之一：数据预处理：caret包应用之一：数据预处理[size=12.7px]
[size=15.3px]在进行数据挖掘时，我们会用到R中的很多扩展包，各自有不同的函数和功能。如果能将它们综合起来应用就会很方便。caret包（Classification and Regression Training）就是为了解决分类和回归问题的数据训练而创建的一个综合工具包。下面的例子围绕数据挖掘的几个核心步骤来说明其应用。
本例涉及到的数据是一个医学实验数据，载入数据之后可以发现其样本数为528，自变量数为342，mdrrDescr为自变量数据框，mdrrClass为因变量。
library(caret)
data(mdrr)本例的样本数据所涉及到的变量非常多，需要对变量进行初步降维。其中一种需要删除的变量是常数自变量，或者是方差极小的自变量，对应的命令是nearZeroVar，可以看到新数据集的自变量减少到了297个。
zerovar=nearZeroVar(mdrrDescr)
newdata1=mdrrDescr[,-zerovar]另一类需要删除的是与其它自变量有很强相关性的变量，对应的命令是findcorrelation。自变量中还有可能存在多重共线性问题，可以用findLinearCombos命令将它们找出来。这样处理后自变量减少为94个。
descrCorr = cor(newdata1)
highCorr = findCorrelation(descrCorr, 0.90)
newdata2 = newdata1[, -highCorr]
comboInfo = findLinearCombos(newdata2)
newdata2=newdata2[, -comboInfo$remove]我们还需要将数据进行标准化并补足缺失值，这时可以用preProcess命令，缺省参数是标准化数据，其高级功能还包括用K近邻和装袋决策树两种方法来预测缺失值。此外它还可以进行cox幂变换和主成分提取。
Process = preProcess(newdata2)
newdata3 = predict(Process, newdata2)最后是用createDataPartition将数据进行划分，分成75%的训练样本和25%检验样本，类似的命令还包括了createResample用来进行简单的自助法抽样，还有createFolds来生成多重交叉检验样本。
inTrain = createDataPartition(mdrrClass, p = 3/4, list = FALSE)
trainx = newdata3[inTrain,]
testx = newdata3[-inTrain,]
trainy = mdrrClass[inTrain]
testy = mdrrClass[-inTrain]在建模前还可以对样本数据进行图形观察，例如对前两个变量绘制箱线图
featurePlot(trainx[,1:2],trainy,plot='box')
caret包应用之二：特征选择[size=12.7px]
[size=15.3px]在进行数据挖掘时，我们并不需要将所有的自变量用来建模，而是从中选择若干最重要的变量，这称为特征选择（feature selection）。一种就是后向选择，即先将所有的变量都包括在模型中，然后计算其效能（如误差、预测精度）和变量重要排序，然后保留最重要的若干变量，再次计算效能，这样反复迭代，找出合适的自变量数目。这种算法的一个缺点在于可能会存在过度拟合，所以需要在此算法外再套上一个样本划分的循环。在caret包中的rfe命令可以完成这项任务。
首先定义几个整数，程序必须测试这些数目的自变量.
subsets = c(20,30,40,50,60,70,80)然后定义控制参数，functions是确定用什么样的模型进行自变量排序，本例选择的模型是随机森林即rfFuncs，可以选择的还有lmFuncs（线性回归），nbFuncs（朴素贝叶斯），treebagFuncs（装袋决策树），caretFuncs（自定义的训练模型）。
method是确定用什么样的抽样方法，本例使用cv即交叉检验, 还有提升boot以及留一交叉检验LOOCV
ctrl= rfeControl(functions = rfFuncs, method = &cv&,verbose = FALSE, returnResamp = &final&)最后使用rfe命令进行特征选择，计算量很大，这得花点时间
Profile = rfe(newdata3, mdrrClass, sizes = subsets, rfeControl = ctrl)观察结果选择50个自变量时，其预测精度最高
print(Profile)Variables Accuracy Kappa AccuracySD KappaSD Selected& && &&&20& &0.5& & 0.50& && && && && &&&30& &0.4& & 0.02& && && && && &&&40& &0.7& & 0.59& && && && && &&&50& &0.6& & 0.44& && &&&*& && &&&60& &0.7& & 0.99& && && && && &&&70& &0.9& & 0.53& && && && && &&&80& &0.9& & 0.73& && && && && &&&94& &0.9& & 0.06用图形也可以观察到同样结果
plot(Profile)
下面的命令则可以返回最终保留的自变量
Profile$optVariables
caret包应用之三：建模与参数优化[size=12.7px]
[size=15.3px]在进行建模时，需对模型的参数进行优化，在caret包中其主要函数命令是train。
首先得到经过特征选择后的样本数据，并划分为训练样本和检验样本
newdata4=newdata3[,Profile$optVariables]
inTrain = createDataPartition(mdrrClass, p = 3/4, list = FALSE)
trainx = newdata4[inTrain,]
testx = newdata4[-inTrain,]
trainy = mdrrClass[inTrain]
testy = mdrrClass[-inTrain]然后定义模型训练参数，method确定多次交叉检验的抽样方法，number确定了划分的重数， repeats确定了反复次数。
fitControl = trainControl(method = &repeatedcv&, number = 10, repeats = 3,returnResamp = &all&)确定参数选择范围，本例建模准备使用gbm算法，相应的参数有如下三项
gbmGrid = expand.grid(.interaction.depth = c(1, 3),.n.trees = c(50, 100, 150, 200, 250, 300),.shrinkage = 0.1)利用train函数进行训练，使用的建模方法为提升决策树方法，
gbmFit1 = train(trainx,trainy,method = &gbm&,trControl = fitControl,tuneGrid = gbmGrid,verbose = FALSE)从结果可以观察到interaction.depth取1，n.trees取150时精度最高
interaction.depth&&n.trees&&Accuracy&&Kappa&&Accuracy SD&&Kappa SD&&1& && && && && && &50& && & 0.822& &&&0.635&&0.0577& && & 0.118& &&&1& && && && && && &100& && &0.824& &&&0.639&&0.0574& && & 0.118& &&&1& && && && && && &150& && &0.826& &&&0.643&&0.0635& && & 0.131& &&&1& && && && && && &200& && &0.824& &&&0.64& &0.0605& && & 0.123& &&&1& && && && && && &250& && &0.816& &&&0.623&&0.0608& && & 0.124& &&&1& && && && && && &300& && &0.824& &&&0.64& &0.0584& && & 0.119& &&&3& && && && && && &50& && & 0.816& &&&0.621&&0.0569& && & 0.117& &&&3& && && && && && &100& && &0.82& && &0.631&&0.0578& && & 0.117& &&&3& && && && && && &150& && &0.815& &&&0.621&&0.0582& && & 0.117& &&&3& && && && && && &200& && &0.82& && &0.63& &0.0618& && & 0.125& &&&3& && && && && && &250& && &0.813& &&&0.617&&0.0632& && & 0.127& &&&3& && && && && && &300& && &0.812& &&&0.615&&0.0622& && & 0.126& &
同样的图形观察
plot(gbmFit1)
caret包应用之四：模型预测与检验[size=12.7px]
[size=15.3px]模型建立好后，我们可以利用predict函数进行预测，例如预测检测样本的前五个
predict(gbmFit1, newdata = testx)[1:5]为了比较不同的模型，还可用装袋决策树建立第二个模型，命名为gbmFit2
gbmFit2= train(trainx, trainy,method = &treebag&,trControl = fitControl)
models = list(gbmFit1, gbmFit2)另一种得到预测结果的方法是使用extractPrediction函数，得到的部分结果如下显示
predValues = extractPrediction(models,testX = testx, testY = testy)
head(predValues)
& &&&obs& &&&pred model dataType&&object1 Active& &Active& &gbm Training Object12 Active& &Active& &gbm Training Object13 Active Inactive& &gbm Training Object14 Active& &Active& &gbm Training Object15 Active& &Active& &gbm Training Object1 6 Active Active gbm Training Object1
从中可提取检验样本的预测结果
testValues = subset(predValues, dataType == &Test&)如果要得到预测概率，则使用extractProb函数
probValues = extractProb(models,testX = testx, testY = testy)
testProbs = subset(probValues, dataType == &Test&)对于分类问题的效能检验，最重要的是观察预测结果的混淆矩阵
Pred1 = subset(testValues, model == &gbm&)
Pred2 = subset(testValues, model == &treebag&)
confusionMatrix(Pred1$pred, Pred1$obs)
confusionMatrix(Pred2$pred, Pred2$obs)结果如下，可见第一个模型在准确率要比第二个模型略好一些
& && && & ReferencePrediction Active Inactive&&Active& && & 65& && & 12&&Inactive& && &9& && & 45& && && && && && && && && && && && && && && && && &&&Accuracy : 0.8397& && && &
ReferencePrediction Active Inactive&&Active& && & 63& && & 12&&Inactive& &&&11& && & 45& && && && && && && && && && && && && && && && && && && &Accuracy : 0.8244& && && &
最后是利用ROCR包来绘制ROC图
prob1 = subset(testProbs, model == &gbm&)
prob2 = subset(testProbs, model == &treebag&)
library(ROCR)
prob1$lable=ifelse(prob1$obs=='Active',yes=1,0)
pred1 = prediction(prob1$Active,prob1$lable)
perf1 = performance(pred1, measure=&tpr&, x.measure=&fpr& )
plot( perf1 )
扫一扫加入本版微信群