刚看完《Introduction to Scientific Programming and Simulation Using R》这本书，
第一部分是R编程入门； 第二部分是数值计算，主要是解方程，求积分和优化；第三部分是概率和统计，主要讲概率、随机变量等概念和参数估计； 第四部分是simulation，主要讲Monte Carlo积分和方差降低。
Monte Carlo方法由冯·诺伊曼于二战时提出，1777年法国人布丰用此思路估计pi值被认为是Monte Carlo的起源，这个方法简单又快速，通过产生随机数，将数值计算问题变成随机变量的某些特征值（比如期望值）。
积分运算，和估计pi值一样，用hit and miss法估计。
hit_miss <- function(fun, lower, upper, n=500) {
# Monte Carlo integration using the hit and miss method
x <- runif(n, lower, upper)
f.value <- sapply(seq(lower, upper, 0.001), fun)
f.min <- min(f.value)
f.max <- max(f.value)
y <- runif(n, f.min, f.max)
hits <- sum(y <= fun(x))
area <- (upper-lower) * f.min + hits/n * (upper-lower) * (f.max-f.min)
return(area)
}
hit and miss方法收敛太慢，效率并不高，通常所说的MC积分是指下面这个方法。
泰勒公式学过微积分都应该知道，可以翻wiki复习一下，https://zh.wikipedia.org/wiki/泰勒公式.
用R简单实现一下：
efv <- function(f, value, variable="x", a=0, eps=0.001) {
#estimate function value using Taylor theorem
assign(eval(variable), a)
fv.old <- eval(f)
k <- 1
repeat {
df <- D(f, variable)
if (df == 0)
break
fv.new <- fv.old + eval(df)*(value-a)^k/factorial(k)
if (fv.new - fv.old < eps)
break
fv.old <- fv.new
f <- df
k <- k+1
}
return (fv.new)
}
翻看了以前写的使用Newton-Raphson Method求一个数的开方，想到其实也可以用中值定理来实现。 中值定理：f(x)是一个连续性的函数，在[u,v]区间内，当c的值位于f(u)和f(v)之间时，至少存在一个点，满足f(x) = c 当f(u)和f(v)一正一负时，那么在[u,v]之间至少有一个根的存在，这个定理本来就是拿来证明根的存在的，但是其实也可以用来求解根。
