R与数学专题 | 第四节: 向量化计算与导数计算
上期内容回顾
上一期我们主要给大家介绍了用R求解线性方程组的解,超定方程组,非线性方程的线性求法。
这一周包括两部分内容:第一部分为向量化计算,主要介绍了向量化计算的思想,并以apply族函数(apply,sapply,lappy,tapply等函数)为例进行了示范;第二部分主要讲解在R中求解导数的方法。
1.向量化计算
向量化计算是一种特殊的并行计算的方式,通常是对不同的数据执行同样的一个或一批指令,或者说把指令应用于一个数组/向量,当我们习惯了用C语言的单步循环的思想来思考问题的时候,要把思维切换到向量化的方式是件比较困难的事件,向量化运算的好处在于避免使用循环,使代码更为简洁、高效和易于理解。R语言以向量为基本运算对象。当输入的对象为向量时,对其中的每个元素分别进行处理,然后以向量的形式输出。R语言中基本上所有的数据运算均能允许向量操作。本文以apply族函数(apply,sapply,lappy,tapply等函数)作一个简单的介绍。
1.1 apply函数
apply函数的处理对象是矩阵或数组,它逐行或逐列的处理数据,其输出的结果将是一个向量或是矩阵。
例如计算矩阵每列的最大值
A=matrix(rpois(16,10),nrow=4,ncol=4,byrow=T) A
## [,1] [,2] [,3] [,4]
## [1,] 12 11 13 6
## [2,] 19 10 11 5
## [3,] 13 16 9 8
## [4,] 4 9 11 3
apply(A,MARGIN=1,max) #MARGIN=1表示按行处理数据,MARGIN=2表示按列
## [1] 13 19 16 11
与按循环的方式程序比对
Max=rep(NA,dim(A)[2])
for(i in 1:dim(A)[1]){
Max[i]=max(A[i,]) }
Max
## [1] 13 19 16 11
从上面两种实现方式,是不是运用apply函数代码更简洁呢。
1.2 lapply函数
lapply的处理对象是向量、列表或其它对象,它将向量中的每个元素作为参数,输入到处理函数中,最后生成结果的格式为列表。 例如要求数据框每列的最大值与最小值:
data <- data.frame(x=rpois(10,10),y=rnorm(10)) data
## x y
## 1 12 0.3536062
## 2 6 -0.3733125
## 3 5 -0.0713025
## 4 9 2.3026481
## 5 6 -0.3563732
## 6 13 0.9635236
## 7 6 0.9434594
## 8 10 0.5080922
## 9 8 0.3922887
## 10 8 -1.5789785
lapply(data,function(data)
c(Max=max(data),Min=min(data)))
## $x
## Max Min
## 13 5
##
## $y
## Max Min
## 2.302648 -1.578979
1.3 sapply函数
sapply在R中使用较频繁,它和lapply是非常相似的,但其输出格式则为矩阵格式。
sapply(data,function(data)
c(Max=max(data),Min=min(data)))
## x y
## Max 13 2.302648
## Min 5 -1.578979
1.4 tapply函数
tapply是专门用来处理分组数据,例如求各科成绩的平均值
data=data.frame(x=rep(c('语文','数学', '外语'),20),y=as.integer(runif(60,60,100)))
head(data)
## x y
## 1 语文 64
## 2 数学 94
## 3 外语 74
## 4 语文 87
## 5 数学 95
## 6 外语 74
然后求解各科成绩的均值:
tapply(data$y,INDEX=data$x,mean)
## 外语 数学 语文
## 78.35 84.75 79.45
通过对这些函数的学习,能够帮助我们更快的处理数据,而不用编辑复杂的函数及使用大量的循环。更多向量化操作函数的使用可以参考《R实战》。
2.导数计算
高等数学是大学学习的一门数学基础课,R语言能实现高等数学的计算。本文主要介绍用R计算导数,可以使用deriv函数直接进行导数的计算。
2.1 初等函数导数的计算
例如:一元多次函数
dy=deriv(y~2*x^2+3*x^3,'x',func=T) dy
## function (x)
## {
## .expr1 <- x^2
## .value <- 2 * .expr1 + 3 * x^3
## .grad <- array(0, c(length(.value), 1L), list(NULL, c("x")))
## .grad[, "x"] <- 2 * (2 * x) + 3 * (3 * .expr1)
## attr(.value, "gradient") <- .grad
## .value
## }
dy(3)
## [1] 99
## attr(,"gradient")
## x
## [1,] 93
正弦函数:
y=sin(x)
dy=deriv(y~sin(x),'x',func=T) dy
## function (x)
## {
## .value <- sin(x)
## .grad <- array(0, c(length(.value), 1L), list(NULL, c("x")))
## .grad[, "x"] <- cos(x)
## attr(.value, "gradient") <- .grad
## .value
## }
dy(pi)
## [1] 1.224647e-16
## attr(,"gradient")
## x
## [1,] -1
2.2 偏导数计算
与一元函数类似,可以这样计算函数的偏导数: 例如:
f(x,y)=x^2+y^2
dxy=deriv(f~x^2+y^2,c('x','y'),func=T) dxy
## function (x, y)
## {
## .value <- x^2 + y^2
## .grad <- array(0, c(length(.value), 2L), list(NULL, c("x",
## "y")))
## .grad[, "x"] <- 2 * x
## .grad[, "y"] <- 2 * y
## attr(.value, "gradient") <- .grad
## .value
## }
dxy(1,1)
## [1] 2
## attr(,"gradient")
## x y
## [1,] 2 2
下面看一个较为复杂的偏导数的例子:

dxy=deriv(f~sin(x*y)+x^y+y^x+exp(x*y),c('x','y'),func=T) dxy
## function (x, y)
## {
## .expr1 <- x * y
## .expr3 <- x^y
## .expr5 <- y^x
## .expr7 <- exp(.expr1)
## .expr9 <- cos(.expr1)
## .value <- sin(.expr1) + .expr3 + .expr5 + .expr7
## .grad <- array(0, c(length(.value), 2L), list(NULL, c("x", "y")))
## .grad[, "x"] <- .expr9 * y + x^(y - 1) * y + .expr5 * log(y) + .expr7 * y
## .grad[, "y"] <- .expr9 * x + .expr3 * log(x) + y^(x - 1) * x + .expr7 * x
## attr(.value, "gradient") <- .grad
## .value
## }
dxy(1,1)
## [1] 5.559753
## attr(,"gradient")
## x y
## [1,] 4.258584 4.258584
计算结果是不是和手算的结果一样啊! 对于手动计算导数时,如果函数比较复杂,那么有时手动计算就会有一定的困难,而用程序计算导数,不会受到公式复杂的影响。
本周小结:
通过本周的学习,大家一定对向量化计算有一定的了解,那么赶快动手起来,把你原来一些程序中涉及循环的程序,都尝试一下运用向量化的改写一下,你一定会发现你程序的性能及可读性有了进一步的提高。
在下一周,将给大家介绍R与计算方法,讲解怎样使自己的程序运行速度更快和更有效率

浙公网安备 33010602011771号